
# set file location here:
dir1 <- "???"

# load an r library, downloading if necesary
library0 <- function(s)
{
	if(s %in% rownames(installed.packages()) == FALSE) {install.packages(s,repos='http://cran.us.r-project.org')}
	require(s,character.only=TRUE)
}

# packages
library0("Formula") # for formulas with multiple parts
library0("fastmatch") # for data lookup
library0("AER") # for IV
library0("numDeriv") # for numerical derivatives
library0("stringr") # for regular expressions
library0("DirectEffects") # for sequential g-estimation

###############
# subroutines #
###############

# vectorized function that determines if a row of a matrix has any missing values
rowmiss <- function(X)
{
	rowSums(is.na(X))>0
}

# get columns while keeping columns names with n x 1 matrix
getrows <- function(A,idx)
{
	if(is.data.frame(A)) {
		B <- as.data.frame(A[idx,])
	} else {
		B <- as.matrix(A[idx,])
	}
	colnames(B) <- colnames(A)
	B
}

# get columns while keeping columns names with n x 1 matrix
getcols <- function(A,idx)
{
	
	if(is.data.frame(A)) {
		B <- as.data.frame(A[,idx])
	} else {
		B <- as.matrix(A[,idx])
	}
	colnames(B) <- colnames(A)[idx]
	B
}

# drop missing observations
drop.na <- function(x) { x[!is.na(x)] }

# get the model matrix, including NAs, from a multipart formula
model.matrix0 <- function(f1,rhs=logical(0),data=NULL)
{
	mf <- model.frame(f1,na.action=na.pass,data=data)
	if(length(rhs) == 0) {
		mm <- model.matrix(f1,data=mf)
	} else {
		mm <- model.matrix(f1,data=mf,rhs=rhs)
	}
	mm
}

# get the model response, including NAs, from a multipart formula
model.response0 <- function(f1,lhs=logical(0),data=NULL)
{
	mf <- model.frame(f1,na.action=na.pass,data=data)
	if(length(lhs) == 0) mr <- model.response(data=mf)
	else                 mr <- model.part(f1,data=mf,lhs=lhs)
	mr
}

# get all the y's and X's from a multipart formula, eliminating NAs and subseting
formula.data <- function(f,subset=logical(0),list1=logical(0),se="",nointercept=FALSE,ind=logical(0),data=NULL)
{
	f1 <- Formula(f)
	len1 <- length(f1)
	if(length(len1) != 2) stop("error in formula.data -- invalid formula")
	y <- list()
	X <- list()
	for(i in 1:len1[1]) y[[i]] <- model.response0(f1,lhs=i,data=data)
	for(i in 1:len1[2]) X[[i]] <- model.matrix0(f1,rhs=i,data=data)
	y0 <- y
	X0 <- X
	n <- nrow(as.matrix(y[[1]]))
	drop <- rep(0,n)
	for(i in 1:len1[1]) drop <- drop | is.na(y[[i]])
	for(i in 1:len1[2]) drop <- drop | rowmiss(X[[i]])

	if(length(subset) == n) {
		drop <- drop | !subset
		drop[is.na(drop)] <- TRUE
	} else if(length(subset) != 0) {
		stop("error in formula data")
	}

	# drop entire individuals
	if(length(ind) > 0) {
		dropind <- ind[drop]
		drop <- drop | !is.na(fmatch(ind,dropind))
	}
	
	if(sum(drop) == n) stop("error in formula.data -- all observations will be dropped")

	for(i in 1:len1[1]) y[[i]] <- y[[i]][!drop]
	for(i in 1:len1[2]) X[[i]] <- getrows(X[[i]],!drop)

	if(nointercept) {
		for(i in 1:len1[2]) {
			if(ncol(X[[i]]) > 1) {
				X[[i]] <- getcols(X[[i]],2:ncol(X[[i]]))
			} else if(ncol(X[[i]]) == 1) {
				X[[i]] <- getcols(X[[i]],integer(0))
			}
		}
		for(i in 1:len1[2]) {
			if(ncol(X0[[i]]) > 1) {
				X0[[i]] <- getcols(X0[[i]],2:ncol(X0[[i]]))
			} else if(ncol(X0[[i]]) == 1) {
				X0[[i]] <- getcols(X0[[i]],integer(0))
			}
		}
	}

	if(length(list1) > 0) for(i in 1:length(list1)) if(length(list1[[i]]) > 0)
	{
		if(is.matrix(list1[[i]])) {
			if(nrow(list1[[i]]) != n) stop("error in formula.data -- element of list1 does not have n rows")
			list1[[i]] <- as.matrix(list1[[i]][!drop,])
		} else if(is.array(list1[[i]])) {
			if(length(dim(list1[[i]])) != 3) stop("error in formula.data -- array is not a 3D array")
			if(dim(list1[[i]])[1] != n) stop("error in formula.data -- element of list1 does not have n rows")
			list1[[i]] <- as.array(list1[[i]][!drop,,])
		} else {
			if(length(list1[[i]]) != n) stop("error in formula.data -- element of list1 is not of length n")
			list1[[i]] <- list1[[i]][!drop]
		}
	}

	ind1 <- ind[!drop]
	for(i in 1:len1[2]) {
		if(ncol(X[[i]]) > 0) {
			keepcol <- colSums(abs(X[[i]]))!=0 # drop columns that are all zero (typically created by factor with subset)
			X[[i]] <- getcols(X[[i]],(1:ncol(X[[i]]))[keepcol])
			X0[[i]] <- getcols(X0[[i]],(1:ncol(X0[[i]]))[keepcol])
		}
	}

	list(y=y,X=X,y0=y0,X0=X0,list1=list1,ind=ind1,samp=(1:length(drop))[!drop])
}

# compute variance-covariance matrix for func1(betahat), where func1 is vector-valued
delta.method <- function(func1,beta,V,...)
{
	c_beta <- jacobian(func1,beta,method="simple",...)
	Vc <- c_beta %*% V %*% t(c_beta)
	list(est=func1(beta,...),se=sqrt(diag(Vc)),V=Vc)
}

# test of a vector of nonlinear restrictions
F.test <- function(func1,N,beta,V,...)
{
	dm <- delta.method(func1,beta,V,...)
	cbeta <- func1(beta,...)
	J <- length(cbeta)
	K <- nrow(V)
	F <- t(cbeta) %*% solve(dm$V) %*% cbeta / J
	pval <- 1 - pf(F,J,N-K)
	list(F=F,df1=J,df2=N-K,pval=pval)
}

# test of a vector of nonlinear restrictions
wald.test <- function(func1,beta,V,...)
{
	dm <- delta.method(func1,beta,V,...)
	cbeta <- func1(beta,...)
	W <- t(cbeta ) %*% solve(dm$V) %*% cbeta
	pval <- 1 - pchisq(W,length(dm$est))
	list(W=W,pval=pval)
}

# is an observations an outlier?
outlier <- function(x,subset=logical(0),cutoff=3)
{
	abs((x - mean0(x)) / sd0(x)) > cutoff
}

# drop outliers
remove.outliers <- function(x,subset=logical(0),cutoff=3)
{
	out <- outlier(x,subset,cutoff)
	if(length(subset) == 0) {
		x[!out]
	} else {
		x[!out&subset]
	}
}

apply0 <- function(func1,x,subset=numeric(0),by=numeric(0),bys=numeric(0),w=logical(0),...)
{
	n <- length(x)
	if(length(subset) != 0 & length(subset) != n) stop("error in apply0 -- incompadible dimensions (subset)")
	if(length(by) != 0 & length(by) != n) stop("error in apply0 -- incompadible dimensions (by)")
	if(length(w) != 0 & length(w) != n) stop("error in apply0 -- incompadible dimensions (w)")

	if(length(subset) == 0) {
		subset1 <- rep(TRUE,n)
	} else {
		subset1 <- subset
	}
	subset1 <- subset1 & !is.na(x)
	if(length(w) != 0) subset1 <- subset1 & !is.na(w)
	subset1[is.na(subset1)] <- FALSE
	if(length(by) != 0) by1 <- subset(by,subset1)
	x1 <- x[subset1]
	if(length(w) != 0) w1 <- w[subset1]
	if(length(x1) == 0) {
		NA
	} else {
		if(length(by) == 0) {
			if(length(w) == 0) {
				func1(x1,...)
			} else {
				func1(x1,w=w1,...)
			}
		} else {
			if(length(bys) == 0) bys <- sort(unique(by1))
		
			# if only one bys is left over, prevent error
			if(length(bys) == 1) {
				if(length(w) == 0) {
					return(func1(x1,...))
				} else {
					return(func1(x1,w=w1,...))
				}
			}
			
			# there has to be a better way than this...
			if(length(w) == 0) {
				dat1 <- cbind(x1,by1)[order(by1),]
			} else {
				dat1 <- cbind(x1,w1,by1)[order(by1),]
			}

			J <- length(unique(by1))
			start <- rep(NA,J)
			start[1] <- 1
			tab1 <- table(by1)

			if(J > 1) {
				for(j in 2:J) start[j] <- start[j-1] + tab1[j-1]
				end <- c(start[2:J]-1,sum(tab1))
			} else {
				end <- sum(tab1)
			}

			ret <- rep(NA,length(bys))
			if(length(w) == 0) {
				for(j in 1:J) {
					k <- fmatch(dat1[start[j],2],bys)
					if(!is.na(k)) ret[k] <- func1(as.numeric(dat1[start[j]:end[j],1]),...)
				}
			} else {
				for(j in 1:J) {
					k <- fmatch(dat1[start[j],2],bys)
					if(!is.na(k)) ret[k] <- func1(as.numeric(dat1[start[j]:end[j],1]),...,w=as.numeric(dat1[start[j]:end[j],2]))
				}
			}
			names(ret) <- bys
			ret
		}
	}
}

mean0 <- function(x,w=numeric(0),subset=numeric(0),by=numeric(0),bys=numeric(0),se=FALSE,ci=FALSE)
{
	n <- length(x)
	if(length(subset) != 0 & length(subset) != n) stop("error in mean0 -- incompadible dimensions (subset)")
	if(length(w) != 0 & length(w) != n) stop("error in mean0 -- incompadible dimensions (w)")
	if(!se & !ci) {
		if(length(w) == 0) {
			apply0(mean,x,subset,by,bys)
		} else {
			if(length(subset) == 0) {
				subset1 <- !is.na(x) & !is.na(w)
			} else {
				subset1 <- !is.na(x) & !is.na(w) & subset
			}
			apply0(mean,w*x,subset1,by,bys) / apply0(mean,w,subset1,by,bys)
		}
	} else {
		ret <- list()
		if(length(w) == 0) {
			ret$est <- apply0(mean,x,subset,by,bys)
			se1 <- ret$se <- (var0(x,subset,by,bys) / n0(x,subset,by,bys))^.5
		} else {
			if(length(subset) == 0) {
				subset1 <- !is.na(x) & !is.na(w)
			} else {
				subset1 <- !is.na(x) & !is.na(w) & subset
			}
			W = apply0(sum,w,subset1,by,bys)
			ret <- list()
			ret$est <- apply0(sum,w*x,subset1,by,bys) / W 
			se1 <- (var0(w*x,subset1,by,bys) * n0(w*x,subset1,by,bys) / W^2)^.5
		}
		if(se) ret$se <- se1
		if(ci) ret$ci <- c(ret$est - 1.96 * se1,ret$est + 1.96 * se1)
		ret
	}
}

quantile0 <- function(x,probs,subset=logical(0),by=numeric(0),bys=numeric(0),w=numeric(0))
{
	n <- length(x)
	if(length(subset) != 0 & length(subset) != n) stop("error in quantile0 -- incompadible dimensions (subset)")
	if(length(w) != 0 & length(w) != n) stop("error in quantile0 -- incompadible dimensions (w)")
	k <- length(probs)
	if(length(by) == 0) {
		m <- 1
	} else {
		if(length(bys) == 0) bys <- sort(unique(by))
		m <- length(bys)
	}
	ret <- matrix(rep(NA,m*k),nrow=m)
	if(length(w) == 0) {
		for(i in 1:k) ret[,i] <- apply0(quantile,x,subset,by,bys,probs=probs[i])
	} else {
		if(length(subset) == 0) subset1 <- rep(TRUE,n)
		subset1 <- subset & !is.na(w)
		for(i in 1:k) ret[,i] <- apply0(wtd.quantile0,x,subset1,by,bys,probs=probs[i],w=w)
	}
	rownames(ret) <- bys
	colnames(ret) <- probs
	ret
}

# sd, removing missing values, and subset option
sd0 <- function(x,subset=numeric(0),by=numeric(0),bys=numeric(0))
{
	if(length(subset) == 0) {
		subset1 <- subset(x,!is.na(x))
		if(length(by) != 0) by1 <- subset(by,!is.na(x))
	} else if(length(subset) != length(x)) {
		# Error
		stop("Incompadible Dimensions")
	} else {
		subset1 <- subset(x,subset&!is.na(x))
		if(length(by) != 0) by1 <- subset(by,subset&!is.na(x))
	}
	if(length(subset1) == 0) {
		NA
	} else {
		if(length(by) == 0) {
			sd(subset1)
		} else {
			if(length(bys) == 0) bys <- sort(unique(by1))
			n <- length(bys)
			ret <- rep(NA,n)
			for(i in 1:n)
			{
				if(length(subset(subset1,by1==bys[i])) > 0) {
					ret[i] <- sd(subset(subset1,by1==bys[i]))
				}
			}
			names(ret) <- bys
			ret
		}
	}
}

# covariance, pairwise
cov0 <- function(x,y,subset=logical(0),by=logical(0),bys=logical(0))
{
	if(length(by) == 0)
	{
		if(length(subset) == 0) {
		cov(x,y,use="pairwise")
		} else {
			x0 <- subset(x,subset)
			y0 <- subset(y,subset)
			cov(x0,y0,use="pairwise")
		}
	} else {
		if(length(bys)==0) bys <- sort(unique(by))
		ret <- rep(NA,length(bys))
		for(i in 1:length(bys)) {
			if(length(subset) == 0) {
				x0 <- subset(x,by==bys[i])
				y0 <- subset(y,by==bys[i])
			} else {
				x0 <- subset(x,subset&by==bys[i])
				y0 <- subset(y,subset&by==bys[i])
			}
			if(length(x0) > 0) {
				ret[i] <- cov(x0,y0,use="pairwise")
			}
		}
		names(ret) <- bys
		ret
	}
}

crossprod0 <- function(x,y)
{
	miss <- !apply(x, 1, function(x)!any(is.na(x))) | !apply(y, 1, function(x)!any(is.na(x)))
	x1 <- x[!miss,]
	y1 <- y[!miss,]
	if(length(x1) == 0) {
		NA
	} else if(length(x1) == ncol(x)) {
		x1 %*% t(y1)
	} else {
		crossprod(x1,y1)
	}
}

# robust standard errors
robust0 <- function(fm)
{
	require(sandwich, quietly = TRUE)
	require(lmtest, quietly = TRUE)
	N <- nrow(estfun(fm))
	M <- N
	K <- length(fm$coef)
	dfc <- (M/(M-1))*((N-1)/(N-K))
	uj  <- apply(estfun(fm),2, function(x) tapply(x, 1:N, sum))
	vcovCL <- dfc*sandwich(fm, meat=crossprod(uj)/N)
	sm <- summary(fm)
	model <- coeftest(fm, vcovCL)
	ret <- list(model=model,N=N,coef=model[1:K,1],theta=model[1:K,1],se=(diag(vcovCL))^.5,V=vcovCL)
	if(!is.null(sm$r.squared)) {
		ret$r.squared <- sm$r.squared
		ret$r.squared.adj = 1 - (1 - ret$r.squared) * (N - K) / (N - 1)
	}
	currclass <- class(fm)
	idx <- fmatch("polr",currclass)
	if(is.na(idx)) idx <- fmatch("tobit",currclass)
	if(is.na(idx)) idx <- fmatch("selection",currclass)
	if(is.na(idx)) idx <- fmatch("glm",currclass)
	if(is.na(idx)) idx <- fmatch("lm",currclass)
	if(!is.na(idx)) {
		currclass <- currclass[idx]
		if(currclass == "lm" | currclass == "glm" | currclass == "selection") {
			ret$log.lik <- logLik(fm)
		} else if(currclass == "polr") {
			ret$log.lik <- logLik.glm(fm)
		} else if(currclass == "tobit") {
			ret$log.lik <- logLik(fm)[1]
		}
	}
	ret
}

# clustered standard errors
cluster0 <- function(fm, cl)
{
	require(sandwich, quietly = TRUE)
	require(lmtest, quietly = TRUE)
	not <- attr(fm$model,"na.action")
	if(!is.null(not)) { cl <- cl[-not] }
	M <- length(unique(cl))
	N <- length(cl)
	K <- fm$rank
	dfc <- (M/(M-1))*((N-1)/(N-K))
	uj <- apply(estfun(fm),2, function(x) tapply(x, cl, sum));
	vcovCL <- dfc*sandwich(fm, meat=crossprod(uj)/N)
	sm <- summary(fm)
	resid <- summary(fm)$residuals
	if(length(resid) > 0) {
		#sigmasqse <- (mean(((M / N) * aggregate((summary(fm)$residuals)^2,by=list(cl),FUN=sum)[2])^2)-(summary(fm)$sigma)^4)^.5 / M^.5
		sigmasqse <- NA
		model <- coeftest(fm, vcovCL)
		ret <- list(model=model,N=N,Clusters=M,sigmasq=summary(fm)$sigma^2,sigmasqse=sigmasqse,coef=sm$coefficients[1:K],theta=sm$coefficients[1:K],se=diag(vcovCL)^.5,V=vcovCL)
	} else {
		model <- coeftest(fm, vcovCL)
		ret <- list(model=model,N=N,Clusters=M,coef=fm$coefficients[1:K],theta=fm$coefficients[1:K],se=diag(vcovCL)^.5,V=vcovCL)
	}
	if(!is.null(sm$r.squared)) ret$r.squared <- sm$r.squared
	currclass <- class(fm)
	idx <- fmatch("polr",currclass)
	if(is.na(idx)) idx <- fmatch("tobit",currclass)
	if(is.na(idx)) idx <- fmatch("selection",currclass)
	if(is.na(idx)) idx <- fmatch("glm",currclass)
	if(is.na(idx)) idx <- fmatch("lm",currclass)
	if(!is.na(idx)) {
		currclass <- currclass[idx]
		if(currclass == "lm" | currclass == "glm" | currclass == "selection") {
			ret$log.lik <- logLik(fm)
		} else if(currclass == "polr") {
			ret$log.lik <- logLik.glm(fm)
		} else if(currclass == "tobit") {
			ret$log.lik <- logLik(fm)[1]
		}
	}
	ret
}

# panel Newey-West standard errors
panelneweywest0 <- function(fm,time,id,lag=logical(0))
{
	not <- attr(fm$model,"na.action")
	if(!is.null(not))
	{
		time <- time[-not]
		id <- id[-not]
	}
	NTot <- length(time)
	uj  <- apply(estfun(fm),2, function(x) tapply(x, 1:NTot, sum))
	K <- ncol(uj)
	T <- max(time) - min(time) + 1
	if(length(lag) == 0) lag <- floor(4 * (T / 100)^(4/9))
	ids <- sort(unique(id))
	times <- min(time):max(time) - min(time) + 1
	N <- length(ids)
	dfc <- NTot / (NTot - K)
	weights <- seq(1, 0, by = -(1/(lag + 1)))
	rval <- matrix(rep(0,K^2),ncol=K)
	
	for(n in 1:N)
	{
		ujcurr <- uj[id==ids[n],]
		timecurr <- time[id==ids[n]] - min(time) + 1
		Tn <- nrow(ujcurr)
		ujcurr2 <- matrix(rep(NA,length(times)*K),ncol=K)
		ujcurr2[timecurr,] <- ujcurr
		dfc <- NTot / (NTot - K)
		weights <- seq(1, 0, by = -(1/(lag + 1)))
		temp <- crossprod0(ujcurr2,ujcurr2)
		if(sum(is.na(temp)) == 0) rval <- rval + 0.5 * weights[1] * temp
  		if(length(weights) > 1)
		{
    			for (ii in 2:length(weights))
			{
				if(sum(is.na(temp)) == 0) rval <- rval + weights[ii] * temp
	    		}
  		}
	}
	rval <- rval + t(rval)
	rval <- dfc  * rval
	rval <- rval / NTot
	vcov <- sandwich(fm, meat=rval)
	model <- coeftest(fm, vcov)
	list(model=model,N=NTot,Clusters=N,T=T,lags=lag,r.squared=summary(fm)$r.squared,coef=model[1:K],theta=model[1:K],se=diag(vcov)^.5,V=vcov)
}


get.ses <- function(m1,se,unitid=logical(0),timeid=logical(0),lag=logical(0))
{
	     if(se == "classic"                     ) stop("error in get.ses --- don't call with classic")
	else if(se == "cluster"                     ) res <- cluster0(m1,unitid)
	else if(se == "nw" | se == "neweywest"      ) res <- neweywest0(m1,timeid)
	else if(se == "pnw" | se == "panelneweywest") res <- panelneweywest0(m1,timeid,unitid,lag)
	else if(se == "robust")                       res <- robust0(m1) # the default option
	else                                          stop("invalid se option")
}

# linear regression
lm0 <- function(f1,se="robust",w=logical(0),subset=logical(0),nointercept=FALSE,unitid=logical(0),timeid=logical(0),lag=logical(0),data=NULL,by=logical(0),bys=logical(0))
{
	if(length(by) > 0) {
		if(length(bys) == 0) bys <- sort(unique(by[subset]))
		lms <- list()
		b <- length(bys)
		for(i in 1:b) {
			if(length(subset) == 0) {
				lms[[i]] <- lm0(f1,se=se,w=w,subset=by==bys[i],nointercept=nointercept,unitid=unitid,timeid=timeid,lag=lag,data=data)
			} else {
				lms[[i]] <- lm0(f1,se=se,w=w,subset=subset&by==bys[i],nointercept=nointercept,unitid=unitid,timeid=timeid,lag=lag,data=data)
			}
		}
		names(lms) <- bys
		return(lms)
	} else {
		# get y and x
		dat <- formula.data(f1,subset,list(unitid=unitid,timeid=timeid,w=w),se=se,nointercept=nointercept,data=data)
		if(length(dat$y) != 1) stop("error in lm0 -- invalid number of repsonses")
		if(length(dat$X) != 1) stop("error in lm0 -- invalid RHS")
		y1 <- dat$y[[1]]
		x1 <- dat$X[[1]]
		w1 <- dat$list1$w
		unitid1 <- dat$list1$unitid
		timeid1 <- dat$list1$timeid

		# estimate the model
		if(length(w1) == 0) {
			lm1 <- lm(y1 ~ 0 + x1)
		} else {
			lm1 <- lm(y1 ~ 0 + x1,weights=w1)
		}

		if(se=="classic") {
			s1 <- summary(lm1)
			res <- list()
			res$coef <- s1$coefficients[,1]
			res$se <- s1$coefficients[,2]
			res$N <- length(y1)
			res$r.squared <- s1$r.squared
			res$ivs <- names(lm1$coefficients)
		} else if(se == "cluster-wild") {
			K <- lm1$rank # number of covariates
			p.val <- lm.cluster.wild(f1,logical(0),x1=x1,y1=y1,unitid=unitid1,forceNull=FALSE)
			res <- get.ses(lm1,se="cluster",unitid=unitid1,lag=lag)
			res$ivs <- rownames(res$model)
			res$p.val <- p.val	
		} else {
			res <- get.ses(lm1,se,timeid=timeid1,unitid=unitid1,lag=lag)
			res$ivs <- rownames(res$model)
		}

		if(length(unitid1) > 0) res$Clusters <- length(unique(unitid1))
		
		# model predictions
		x <- dat$X0[[1]]
		pred <- as.vector(x %*% lm1$coef)
		XpX <- t(x1) %*% x1
		res$pred <- pred
		res$res <- dat$y0[[1]][[1]] - pred
		res$XpX <- XpX

		# hack to get correct model fit when intercept is included as a variable
		if(length(w1) == 0) {
			lm2 <- lm(y1 ~ x1)
		} else {
			lm2 <- lm(y1 ~ x1,weights=w1)
		}

		res$theta <- res$coef
		res$r.squared <- summary(lm2)$r.squared
		res$sigma <- summary(lm2)$sigma
		res$AIC <- AIC(lm2)
		res$BIC <- BIC(lm2)
		fac <- (res$N-1)/(res$N-ncol(x))
		res$r.squared.adj <- fac * res$r.squared + 1 - fac
		res$ivs <- substr(res$ivs,3,nchar(res$ivs))
		res$ivs[res$ivs=="tercept)"] <- "Intercept"
		if(se=="classic") {
			res$V <- res$sigma^2 * solve(XpX)
			res$predse <- rep(NA,nrow(x))
			for(n in 1:nrow(x)) res$predse[n] <- res$sigma * (1 + t(x[n,]) %*% solve(XpX) %*% x[n,])^.5
			tcrit <- qt(.975,res$N-ncol(x))
			res$predlow <- res$pred - tcrit * res$predse
			res$predhigh <- res$pred + tcrit * res$predse
		}
		res
	}
}

# IV regression, with subset and cluster options
ivreg0 <- function(f1,se="robust",subset=logical(0),unitid=logical(0),timeid=logical(0),liml=FALSE,lag=logical(0),w=logical(0),data=NULL)
{
	# get y, X, and Z
	dat <- formula.data(f1,subset,list(unitid=unitid,timeid=timeid,w=w),se=se,nointercept=TRUE,data=data)
	if(length(dat$y) != 1) stop("error in ivreg0 -- invalid number of repsonses")
	if(length(dat$X) != 2) stop("error in ivreg0 -- invalid RHS")
	y <- dat$y[[1]]
	X <- as.matrix(dat$X[[1]])
	Z <- as.matrix(dat$X[[2]])
	if(length(w) != 0) W <- dat$list1$w
	unitid1 <- dat$list1$unitid
	timeid1 <- dat$list1$timeid
	XX <- getcols(X,is.na(fmatch(colnames(X),colnames(Z)))) # endogenous
	ZZ <- getcols(Z,is.na(fmatch(colnames(Z),colnames(X)))) # instr
	XZ <- getcols(X,!is.na(fmatch(colnames(X),colnames(Z)))) # exogenous
	L <- ncol(XX) # number of endo.

	# estimate the model
	res <- list()
	if(!liml) {
		# 2SLS
		if(length(w) == 0) {
			iv1 <- ivreg(y ~ X,~Z)
		} else {
			iv1 <- ivreg(y ~ X,~Z,weights=W)
		}
		res <- get.ses(iv1,se,timeid=timeid1,unitid=unitid1,lag=lag)
		if(length(unitid1) > 0) {
			res$Clusters <- length(unique(unitid1))
		}
	} else {
		if(length(w) != 0) stop("error in ivreg0 -- liml doesn't work with weights yet")
		# LIML
		Z1 <- cbind(1,Z)
		X1 <- cbind(1,X)
		XZ1 <- cbind(1,XZ)
		N <- nrow(Z1)
		R <- cbind(y,XX)
		M  <- diag(N) - Z1 %*% solve((t(Z1) %*% Z1)) %*% t(Z1)	
		Mi <- diag(N) - XZ1 %*% solve((t(XZ1) %*% XZ1)) %*% t(XZ1)
		S <- t(R) %*% Mi %*% R %*% solve(t(R) %*% M %*% R)
		lambda <- min(eigen(S)$values)
		res$N <- N
		bread <- solve(t(X1) %*% (diag(N) - lambda * M) %*% X1)
		res$theta <- as.vector(bread %*% t(X1) %*% (diag(N) - lambda * M) %*% y)
	
		eps <- y - X1 %*% res$theta
		v <- as.vector(eps) * ((diag(N) - lambda * M) %*% X1)
		if(se == "robust") {
			meat <- var(v)
			res$V <- N * bread %*% meat %*% bread
		} else {
			I <- length(unique(unitid1))
			v2 <- matrix(rep(NA,I*ncol(v)),nrow=I)
			for(j in 1:ncol(v)) {
				v2[,j] <- aggregate(v[,j],by=list(unitid1),FUN=sum)[[2]]
			}
			meat <- var(v2)
			res$V <- I * bread %*% meat %*% bread
			res$Clusters <- I
		}
		res$se <- diag(res$V)^.5
	}
	res$ivs <- c("(Intercept)",colnames(X)) # rename for consistency
	
	# weak instrument check:
	for(l in 1:L) {
		if(ncol(XZ) > 0) {
			if(length(w) > 0) {
				lm1 <- lm0(XX[,l] ~ ZZ + XZ,se="classic",w=W)
			} else {
				lm1 <- lm0(XX[,l] ~ ZZ + XZ,se="classic")
			}
		} else {
			if(length(w) > 0) {
				lm1 <- lm0(XX[,l] ~ ZZ,se="classic",w=W)
			} else {
				lm1 <- lm0(XX[,l] ~ ZZ,se="classic")
			}
		}
		res[[paste("FirstStageF_",l,sep="")]] <- F.test(function(x) x[fmatch(paste("ZZ",colnames(ZZ),sep=""),lm1$ivs)],lm1$N,lm1$coef,lm1$V)$F
	}
	res$r.squared <- NULL
	res
}

# simulate p-value for a weighted sum of chi-sqaured rv's
wchi_pval <- function(v1,zeta,S)
{
	m <- length(zeta)
	sum((matrix(rnorm(S*m)^2,nrow=S) %*% zeta) > v1) / S
}

vuongtest.lm.func1 <- function(theta,y,X)
{
	beta <- theta[1:(length(theta)-1)]
	sigma <- theta[length(theta)]
	log(exp(-0.5 * ((y - X  %*% beta) / sigma)^2) / (sqrt(2 * pi) * sigma))
}

vuongtest.lm.func2 <- function(theta,y,X)
{
	beta <- theta[1:(length(theta)-1)]
	sigma <- theta[length(theta)]
	-sum(log(exp(-0.5 * ((y - X %*% beta) /sigma)^2) / (sqrt(2 * pi) * sigma)))
}

vuongtest.lm <- function(f,S=1000,alpha=0.05,subset=logical(0),trans=NULL,data=NULL,se="robust",timeid=logical(0))
{
	dat <- formula.data(f,subset,data=data)
	if(length(dat$y) != 1) stop("error in vuongtest.lm -- invalid number of repsonses")
	if(length(dat$X) == 1) {
		if(is.null(trans)) {
			stop("error in vuongtest.lm -- invalid RHS when trans is missing")
		}
	} else if(length(dat$X) != 2) {
		stop("error in vuongtest.lm -- invalid RHS")
	}

	n <- length(dat$y[[1]])
	if(!is.null(trans)) {
		if(!is.function(trans)) {
			stop("error in vuongtest.lm -- trans is not a function")
		}
	}
	y1 <- dat$y[[1]]
	if(is.null(trans)) {
		y2 <- y1
	} else {
		y2 <- trans(y1)
	}
	if(sum(is.na(y2))) {
		stop("error in vuongtest.lm -- trans introduce missing values")
	}
	X1 <- dat$X[[1]]
	if(length(dat$X) == 1) {
		X2 <- X1
	} else {
		X2 <- dat$X[[2]]
	}

	# estimates
	lm1 <- lm0(y1 ~ 0 + X1)
	lm2 <- lm0(y2 ~ 0 + X2)
	beta <- lm1$coef
	gamma <- lm2$coef
	beta <- c(beta,mean0((y1 - X1 %*% beta)^2)^.5)
	gamma <- c(gamma,mean0((y2 - X2 %*% gamma)^2)^.5)

	# compute LR-statistic
	logf <- vuongtest.lm.func1(beta,y1,X1)
	logg <- vuongtest.lm.func1(gamma,y2,X2)
	LR <- sum(logf - logg)

	# derivatives
	dlogffunc <- numDeriv::jacobian(vuongtest.lm.func1,x=beta,y=y1,X=X1)
	dloggfunc <- numDeriv::jacobian(vuongtest.lm.func1,x=gamma,y=y2,X=X2)
	hessf <- numDeriv::hessian(vuongtest.lm.func2,x=beta,y=y1,X=X1)
	hessg <- numDeriv::hessian(vuongtest.lm.func2,x=gamma,y=y2,X=X2)

	# A's stay the same with newey west
	Af <- -hessf / n 
	Ag <- -hessg / n
	AfInv <- solve(Af)
	AgInv <- solve(Ag)

	# B's and omega2 depend on type of covariance estimator
	# W does not reduce in the same way with nw becasue information equality does not apply
	if(se == "robust") {
		Bf <- t(dlogffunc) %*% dlogffunc / n
		Bg <- t(dloggfunc) %*% dloggfunc / n
		Bfg <- t(dlogffunc) %*% dloggfunc / n
		omega2 <- sum((logf - logg)^2) / n - (LR / n)^2
		W11 <- -Bf %*% AfInv
		W12 <- -Bfg %*% AgInv
		W21 <- t(Bfg) %*% AfInv
		W22 <- Bg %*% AgInv
		W <- rbind(cbind(W11,W12),cbind(W21,W22))
	} else if(se == "nw" | se == "neweywest") {
		K1 <- ncol(dlogffunc)
		K2 <- ncol(dloggfunc)
		B <- neweywest.var(cbind(dlogffunc,dloggfunc))
		omega2 <- neweywest.var(logf - logg)
		A <- bdiag(AfInv,AgInv)
		W <- B %*% A
	} else {
		print0("error in vuongtest.lm -- invalid se option")
	}
	zeta <- (eigen(W)$values)^2

	v1 <- n * omega2 # first Vuong test statistic -- is parameter in overlapping region?
	v2 <- LR / (omega2^.5 * n^.5) # second Vuong test statastic -- do the models fit equally well?
	pval1 <- wchi_pval(v1,zeta,S)
	pval2 <- 2 * pnorm(-abs(v2))

	if(pval1 >= alpha)
	{
		message <- "accept null of parameter in overlapping region"
	} else if(pval2 >= alpha) {
		message <- "accept null of equal fit"
	} else if(v2 > 0) {
		message <- "model 1 is better"
	} else {
		message <- "model 2 is better"
	}
	list(message=message,pval1=pval1,pval2=pval2)
}

vuongtest.lm.multi <- function(f,S=1000,alpha=0.05,subset=logical(0),data=NULL,se="robust",timeid=logical(0),digk=3)
{
	dat <- formula.data(f,subset,data=data)
	if(length(dat$y) != 1) stop("error in vuongtest.lm -- invalid number of repsonses")
	n <- length(dat$y[[1]])
	J <- length(dat$X)

	# initial values
	lms <- list()
	betas <- list()
	r.squared <- rep(NA,J)
	dlogffuncs <- list()	
	hesss <- list()
	Bs <- list()
	AInvs <- list()
	logfs <- list()
	for(j in 1:J)
	{
		lms[[j]] <- lm(dat$y[[1]] ~ dat$X[[j]][,2:ncol(dat$X[[j]])])
		r.squared[j] <- dig3(summary(lms[[j]])$r.squared)
		betas[[j]] <- lms[[j]]$coef
		betas[[j]] <- c(betas[[j]],mean0((dat$y[[1]] - dat$X[[j]] %*% betas[[j]])^2)^.5)
		dlogffuncs[[j]] <- numDeriv::jacobian(vuongtest.lm.func1,x=betas[[j]],y=dat$y[[1]],X=dat$X[[j]])
		hesss[[j]] <- numDeriv::hessian(vuongtest.lm.func2,x=betas[[j]],y=dat$y[[1]],X=dat$X[[j]])	
		if(se == "robust") {
			Bs[[j]] <- t(dlogffuncs[[j]]) %*% dlogffuncs[[j]] / n
		} else if(se == "nw" | se == "neweywest") {
			# do nothing
		} else {
			print0("error in vuongtest.lm.multi -- invalid se option")
		}
		AInvs[[j]] <- solve(-hesss[[j]] / n)
		logfs[[j]] <- vuongtest.lm.func1(betas[[j]],dat$y[[1]],dat$X[[j]])
	}

	message <- matrix(rep(NA,J*J),nrow=J)
	pval1 <- matrix(rep(NA,J*J),nrow=J)
	pval2 <- matrix(rep(NA,J*J),nrow=J)

	for(j in 1:(J-1)) for(k in (j+1):J)
	{
		LR <- sum(logfs[[j]] - logfs[[k]])
		if(se == "robust") {
			Bfg <- t(dlogffuncs[[j]]) %*% dlogffuncs[[k]] / n
			W11 <- -Bs[[j]] %*% AInvs[[j]]
			W12 <- -Bfg %*% AInvs[[k]]
			W21 <- t(Bfg) %*% AInvs[[j]]
			W22 <- Bs[[k]] %*% AInvs[[k]]
			W <- rbind(cbind(W11,W12),cbind(W21,W22))
			omega2 <- sum((logfs[[j]] - logfs[[k]])^2) / n - (LR / n)^2
		} else if(se == "nw" | se == "neweywest") {
			K1 <- ncol(dlogffuncs[[j]])
			K2 <- ncol(dlogffuncs[[k]])
			B <- neweywest.var(cbind(dlogffuncs[[j]],dlogffuncs[[k]]))
			A <- bdiag(AInvs[[j]],AInvs[[k]])
			W <- B %*% A
			omega2 <- neweywest.var(logfs[[j]] - logfs[[k]])
		} else {
			print0("error in vuongtest.lm.multi -- invalid se option")
		}

		zeta = (eigen(W)$values)^2
		v1 <- n * omega2 # first Vuong test statistic
		v2 <- LR / (omega2^.5 * n^.5)
		pval1[j,k] <- pval.as.pval(wchi_pval(v1,zeta,S),k=digk)
		pval2[j,k] <- pval.as.pval(2 * pnorm(-abs(v2)),k=digk)
		if(pval1[j,k] >= alpha) {
			message[j,k] <- "accept null of parameter in overlapping region"
		} else if(pval2[j,k] >= alpha) {
			message[j,k] <- "accept null of equal fit"
		} else if(v2 > 0) {
			message[j,k] <- "row model is better"
		} else {
			message[j,k] <- "column model is better"
		}

		pval1[k,j] <- pval1[j,k]
		pval2[k,j] <- pval2[j,k]
		if(message[j,k] == "row model is better") {
			message[k,j] <- "column model is better"
		} else if(message[j,k] == "column model is better") {
			message[k,j] <- "row model is better"
		} else {
			message[k,j] <- message[j,k]
		}
	}

	rownames(message) <- paste("model ",1:J,sep="")
	colnames(message) <- paste("model ",1:J,sep="")
	rownames(pval1) <- paste("model ",1:J,sep="")
	colnames(pval1) <- paste("model ",1:J,sep="")
	rownames(pval2) <- paste("model ",1:J,sep="")
	colnames(pval2) <- paste("model ",1:J,sep="")
	list(r.squared=r.squared,message=message,pval1=pval1,pval2=pval2)
}

dig <- function(x,k=3)
{
	y <- format(round(x,k),nsmall=k,trim=T,scientific=F)
	for(i in 1:length(y)) {
		if(is.na(x[i])) y[i] <- ""
	}
	y
}

dig0 <- function(x)
{
	dig(x,0)
}

dig3 <- function(x)
{
	dig(x,3)
}

as.coef <- function(coef,se,k=3,p.val=logical(0),sigbase=logical(0))
{
	if(length(p.val) == 0) {
		if(is.matrix(coef)) {
			res <- coef
			if(length(sigbase) == 0) {
				for(i in 1:ncol(coef)) res[,i] <- paste(dig(coef[,i],k=k),stars(coef[,i],se[,i]),sep="")
			} else {
				for(i in 1:ncol(coef)) res[,i] <- paste(dig(coef[,i],k=k),stars(coef[,i],se[,i],sigbase=sigbase[i]),sep="")
			}
			res
		} else {
			paste(dig(coef,k=k),stars(coef,se,sigbase=sigbase),sep="")
		}
	} else {
		if(is.matrix(coef)) {
			res <- coef
			if(length(sigbase) == 1) {
				for(i in 1:ncol(coef)) res[,i] <- paste(dig(coef[,i],k=k),stars.pval(coef[,i],p.val[,i],sigbase=sigbase),sep="")
			} else {
				for(i in 1:ncol(coef)) res[,i] <- paste(dig(coef[,i],k=k),stars.pval(coef[,i],p.val[,i],sigbase=sigbase[i]),sep="")
			}
			res
		} else {
			paste(dig(coef,k=k),stars.pval(p.val,sigbase=sigbase),sep="")
		}
	}
}

as.se <- function(se,k=3)
{
	if(is.matrix(se)) {
		res <- se
		for(i in 1:ncol(se)) res[,i] <- str_replace_all(paste("(",dig(se[,i],k),")",sep=""),"\\(\\)","")
		res
	} else {
		str_replace_all(paste("(",dig(se,k),")",sep=""),"\\(\\)","")
	}
}

pval <- function(est,se,sigbase=logical(0))
{
	if(length(sigbase) == 0) {
		2*pnorm(-abs(est/se))
	} else if(length(sigbase) == 1) {
		2*pnorm(-abs((est-sigbase)/se))
	} else {
		stop("error in pval -- length(sigbase) != 0 or 1")
	}
}

stars.scalar <- function(coef,se,sigbase=logical(0),sigalt="^")
{
	if(is.na(coef)) {
		""
	} else if(is.na(se)) {
		"?"
	} else {
		pval1 <- pval(coef,se,sigbase)
		if(is.nan(pval1)) {
			"?"
		} else {
			check1 <- length(sigbase) == 0
			if(!check1) check1 <- sigbase == 0
			if(check1) {
				if(pval1 <= 0.001) {
					"***"
				} else if(pval1 <= 0.01) {
					"**"
				} else if(pval1 <= 0.05) {
					"*"
				} else if(pval1 <= 0.1) {
					"+"
				} else {
					""
				}
			} else {
				if(pval1 <= 0.05) {
					sigalt
				} else {
					""
				}
			}
		}
	}
}

stars <- function(coef,se,sigbase=logical(0))
{
	if(is.matrix(coef))
	{
		m <- nrow(coef)
		n <- ncol(coef)
		res <- coef
		if(length(sigbase) == 0) {
			for(i in 1:m) for(j in 1:n) res[i,j] <- stars.scalar(coef[i,j],se[i,j])
		} else {
			for(i in 1:m) for(j in 1:n) res[i,j] <- stars.scalar(coef[i,j],se[i,j],sigbase[j])
		}
	} else {
		m <- length(coef)
		res <- coef
		if(length(sigbase) == 0) {
			for(i in 1:m) res[i] <- stars.scalar(coef[i],se[i])
		} else if(length(sigbase) == 1) {
			for(i in 1:m) res[i] <- stars.scalar(coef[i],se[i],sigbase=sigbase)
		} else {
			for(i in 1:m) res[i] <- stars.scalar(coef[i],se[i],sigbase=sigbase[i])
		}
	}
	res
}

stars.pval <- function(pval1,sigbase=logical(0))
{
	n <- length(pval1)
	ret <- rep(NA,n)
	check1 <- length(sigbase) == 0
	if(!check1) check1 <- sigbase == 0
	if(check1) {
		for(i in 1:n) {
			if(pval1[i] <= 0.001) {
				ret[i] <- "***"
			} else if(pval1[i] <= 0.01) {
				ret[i] <- "**"
			} else if(pval1[i] <= 0.05) {
				ret[i] <- "*"
			} else if(pval1[i] <= 0.1) {
				ret[i] <- "+"
			} else {
				ret[i] <- ""
			}
		}
	} else {
		for(i in 1:n) {
			if(pval1[i] <= 0.001) {
				ret[i] <- "^^^"
			} else if(pval1[i] <= 0.01) {
				ret[i] <- "^^"
			} else if(pval1[i] <= 0.05) {
				ret[i] <- "^"
			} else {
				ret[i] <- ""
			}
		}
	}
	ret
}

pval.as.pval <- function(pval,k=3)
{
	paste(dig(pval,k),stars.pval(pval),sep="")
}

panel.lag <- function(x,time,id,lags=1,vec=F)
{
	n <- length(x)
	ids <- sort(unique(id))
	m <- length(ids)
	if(vec) {
		if(length(lags) != n) stop("error in panel.lag -- incompadible dimensions")
		xlag <- matrix(rep(NA,n*1),ncol=1)
		for(i in 1:n) {
			k <- fmatch(time[i]-lags[i],time)
			if(!is.na(k)) xlag[i,1] <- x[k]
		}
	} else {
		p <- length(lags)
		xlag <- matrix(rep(NA,n*p),ncol=p)
		xs <- list()
		times <- list()
		for(i in 1:m)
		{
			xs[[i]] <- x[id==ids[i]]
			times[[i]] <- time[id==ids[i]]
		}
		for(i in 1:n)
		{
			j <- fmatch(id[i],ids)
			k <- fmatch(time[i]-lags,times[[j]])
			for(l in 1:length(k))
			if(!is.na(k[l])) {
				xlag[i,l] <- xs[[j]][k[l]]
			}
		}
	}
	xlag
}

output.models <- function(ms,stats=c("N","T","Clusters","r.squared"),rename=logical(0))
{
	J <- length(ms)

	# rename ivs and stats
	if(length(rename) > 0) {
		for(j in 1:J) {
			for(k in 1:length(rename)) {
				# ivs
				idx <- fmatch(names(rename)[k],ms[[j]]$ivs)
				if(!is.na(idx)) ms[[j]]$ivs[idx] <- rename[[k]]

				# stats
				idx <- fmatch(names(rename)[k],names(ms[[j]]))
				if(!is.na(idx)) names(ms[[j]])[idx] <- rename[[k]]

				if(length(stats) > 0) {
					idx <- fmatch(names(rename)[k],stats)
					if(!is.na(idx)) stats[idx] <- rename[[k]]
				}
			}
		}
	}

	ivs <- ms[[1]]$ivs
	if(J > 1) for(j in 2:J) ivs <- c(ivs,ms[[j]]$ivs)
	ivs <- unique(ivs)

	# reorder ivs
	if(length(rename) > 0) {
		newnames <- as.vector(rename)
		newnames <- intersect(newnames,ivs)
		oldnames <- setdiff(ivs,newnames)
		ivs <- c(newnames,oldnames)
	}

	K <- length(ivs)
	coef <- matrix(rep(NA,J*K),ncol=J)
	se <- matrix(rep(NA,J*K),ncol=J)

	if(length(stats) > 0) {
		stats1 <- logical(0)
		for(j in 1:J) {
			idx <- drop.na(fmatch(stats,names(ms[[j]])))
			stats1 <- c(stats1,names(ms[[j]])[idx])
		}
		stats1 <- unique(stats1)
		stats1 <- stats[!is.na(fmatch(stats,stats1))]
		extra <- matrix(rep(NA,J*length(stats1)),ncol=J) 
	}

	for(j in 1:J)
	{
		for(k in 1:length(ms[[j]]$ivs))
		{
			coef[fmatch(ms[[j]]$ivs[k],ivs),j] <- ms[[j]]$theta[k]
			se[fmatch(ms[[j]]$ivs[k],ivs),j] <- ms[[j]]$se[k]
		}

		if(length(stats > 0)) {
			for(i in 1:length(stats1)) {
				idx <- fmatch(stats1[i],names(ms[[j]]))
				if(!is.na(idx)) extra[i,j] <- ms[[j]][[idx]]
			}
		}
	}

	rownames(coef) <- ivs
	colnames(coef) <- names(ms)
	rownames(se) <- ivs
	colnames(se) <- names(ms)
	if(length(stats) > 0) {
		rownames(extra) <- stats1
		list(coef=coef,se=se,stats=extra)
	} else {
		list(coef=coef,se=se)
	}
}

pad0 <- function(x)
{
	K <- length(x)
	x[is.na(x)] <- ""
	len <- max(nchar(x))
	for(k in 1:K) if(nchar(x[k]) < len) x[k] <- paste(x[k],paste(rep(" ",len-nchar(x[k])),collapse=""),sep="")
	x
}

print.output <- function(o,kcoef=3,kstat=3)
{
	J <- ncol(o$coef)
	ivs <- rownames(o$coef)
	stats <- rownames(o$stats)	
	v1 <- c("","Variable",ivs,stats)
	K <- nrow(o$stats)
	dig_stats <- dig(o$stats,k=kstat)
	for(i in 1:nrow(dig_stats)) {
		#if(rownames(dig_stats)[i] == "N" | rownames(dig_stats)[i] == "Clusters" | rownames(dig_stats)[i] == "C" | rownames(dig_stats)[i] == "Observations") {
			dig_stats[i,] <- as.character(as.numeric(dig_stats[i,]))
		#}
	}
	v2 <- rbind(colnames(o$coef),rep("Coef.",J),as.coef(o$coef,o$se,k=kcoef),dig_stats)
	v3 <- rbind(rep("",J),rep("(s.e.)",J),as.se(o$se,k=kcoef),matrix(rep(NA,J*K),nrow=K))
	v1 <- pad0(v1)
	for(j in 1:J) {
		v2[,j] <- pad0(v2[,j])
		v3[,j] <- pad0(v3[,j])
		v2[,j] <- paste(v2[,j],v3[,j],sep=" ")
	}
	for(i in 1:length(v1)) print0(paste(c(v1[i],v2[i,]),collapse=" "))
}

alternate <- function(x,y)
{
	if(is.matrix(x)) {
		m <- nrow(x)
		n <- ncol(x)
		res <- matrix(rep(NA,length(x)*2),ncol=n)
		for(i in 1:n) {
			res[,i] <- as.vector(t(matrix(c(x[,i],y[,i]),nrow=m)))
		}
		res
	} else {
		n <- length(x)
		if(length(y) != n) {
			error("Error with alternate --- incompadible dimensions")
		}
		as.vector(t(matrix(c(x,y),nrow=n)))
	}
}

print0 <- function(s1,s2=logical(0))
{
	if(length(s2) == 0)
	{
		print(s1)
		flush.console()
	} else {
		print(paste(s1,": ",s2,sep=""))
		flush.console()
	}
}

points0 <- function(x,y,subset=logical(0),type='p',col='black',addreg=FALSE,pch=1,lty=1,regcol='black',lwd=NULL,xse=logical(0),yse=logical(0),sort=T,smoothonly=F,bw=logical(0),addsmooth=F,ciopt=1)
{
	# auto sort by x
	if(sort) {
		#if(type=='l' & length(labels) == 0 & length(xse)==0 & length(yse)==0) {
		if(type=='l' & length(xse)==0 & length(yse)==0) {
			if(length(subset) > 0) {
				temp <- cbind(x,y,subset)
				temp <- temp[order(x),]
				x <- temp[,1]
				y <- temp[,2]
				subset <- as.logical(temp[,3])
			} else {
				temp <- cbind(x,y)
				temp <- temp[order(x),]
				x <- temp[,1]
				y <- temp[,2]
			}
		}
	}

	if(length(subset) == 0) {
		if(smoothonly) {
			points(ksmooth0(x,y,bandwidth=bw),type=type,col=col,pch=pch,lty=lty,lwd=lwd)
		} else {
			if(ciopt == 1 & length(yse) > 0) {
				points.ci(x,y,yse,col=col,lty=lty,lwd=lwd)
			} else {
				points(x,y,type=type,col=col,pch=pch,lty=lty,lwd=lwd)
			}
			if(ciopt == 2) {
				if(length(xse) > 0) segments(x-1.96*xse,y,x1=x+1.96*xse,y1=y)
				if(length(yse) > 0) segments(x,y-1.96*yse,x,y1=y+1.96*yse)
			}
			if(addreg) { abline(lm(y ~ x),col=regcol) }
			if(addsmooth) { points(ksmooth0(x,y,bandwidth=bw),type='l',col='magenta') }
		}
	} else {
		x0 <- subset(x,subset)
		y0 <- subset(y,subset)
		yse0 <- subset(yse,subset)
		if(smoothonly) {
			points(ksmooth0(x0,y0,bandwidth=bw),type=type,col=col,pch=pch,lty=lty,lwd=lwd)
		} else {
			if(ciopt == 1 & length(yse) > 0) {
				points.ci(x0,y0,yse0,col=col,lty=lty,lwd=lwd)
			} else {
				points(x0,y0,type=type,col=col,pch=pch,lty=lty,lwd=lwd)
			}
			if(ciopt == 2) {
				if(length(xse) > 0) segments(x0-1.96*xse,y0,x1=x0+1.96*xse,y1=y0)
				if(length(yse) > 0) segments(x0,y0-1.96*yse,x0,y1=y0+1.96*yse)
			}
			if(addreg) { abline(lm(y0 ~ x0),col=regcol) }
			if(addsmooth) { points(ksmooth0(x0,y0,bandwidth=bw),type='l',col='magenta') }
		}
	}
}

latex_doc_start <- function(file1,packages=logical(0))
{
	cat("\\documentclass[11pt]{article}",file=file1,sep="\n")
	s <- "\\usepackage{float,amsmath,amssymb,amsthm,fullpage,setspace,natbib,enumerate,url,multirow,lscape"
	if(length(packages) > 0) for(i in 1:length(packages)) s <- paste(s,",",packages[i],sep="")
	s <- paste(s,"}",sep="")
	cat(s,file=file1,sep="\n",append=TRUE)
	#cat("\\usepackage[dvips]{graphicx,xcolor}",file=file1,sep="\n",append=TRUE)
	cat("\\onehalfspacing",file=file1,sep="\n",append=TRUE)
	cat("\\begin{document}",file=file1,sep="\n\n",append=TRUE)
}

latex_doc_end <- function(file1)
{
	cat("\\end{document}",file=file1,sep="\n",append=TRUE)
}

latex_table_start <- function(file1,k)
{
	cat("\\begin{table}[htp]\\centering\\scriptsize",file=file1,sep="\n",append=TRUE)
	cat("\\begin{singlespace}",file=file1,sep="\n",append=TRUE)
	cat("\\begin{tabular}{ l ",paste(rep("c",k),collapse=" ")," } \\hline \\hline",file=file1,sep="\n",append=TRUE)
}

latex_table_end <- function(file1,caption="",label="")
{
	cat("\\hline \\hline",file=file1,sep="\n",append=TRUE)
	cat("\\end{tabular}",file=file1,sep="\n",append=TRUE)
	cat("\\end{singlespace}",file=file1,sep="\n",append=TRUE)
	if(caption != "") cat(paste("\\caption{\\footnotesize{",caption,"}}",sep=""),file=file1,sep="\n",append=TRUE)
	if(label != "") cat(paste("\\label{",label,"}",sep=""),file=file1,sep="\n",append=TRUE)
	cat("\\end{table}",file=file1,sep="\n",append=TRUE)
}

latex_add_line <- function(file1,s)
{
	cat(s,file=file1,sep="\n",append=TRUE)
}

latex_table_entry <- function(file1,name,coef,bold=FALSE,italic=FALSE)
{
	if(bold)
	{
		name <- paste("\\textbf{",name,"}",sep="")
		coef <- paste("\\textbf{",coef,"}",sep="")
	}
	if(italic)
	{
		name <- paste("\\textit{",name,"}",sep="")
		coef <- paste("\\textit{",coef,"}",sep="")
	}
	s <- name
	n <- length(coef)
	for(i in 1:n) if(is.na(coef[i])) { coef[i] <- "" }
	if(n > 1) for(i in 1:(n-1)) { s <- paste(s," & ",coef[i],sep="") }
	s <- paste(s," & ",coef[n]," \\\\ \n",sep="")
	s
	cat(s,file=file1,sep="\n",append=TRUE)
}

latex_table_entries <- function(file1,names,coefs,bold=FALSE,italic=FALSE)
{
	m <- length(names)
	for(j in 1:m) latex_table_entry(file1,names[j],coefs[j,],bold)
}

latex_model <- function(file1,names,o1,kcoef=3,kstat=3,eqlabels=TRUE,eqlabelsskip=TRUE,stats=c("N","T","Clusters","r.squared"),sigbase=logical(0),modellab="")
{
	m <- length(names)
	n <- ncol(o1$coef)
	if(length(stats) > 0) {
		stats1 <- logical(0)
		stats1 <- drop.na(fmatch(stats,rownames(o1$stats)))
	}
	
	if(eqlabels) {
		if(modellab == "") {
			latex_table_entry(file1,"",colnames(o1$coef))
		} else {
			latex_table_entry(file1,paste("\\textbf{",modellab,":}",sep=""),colnames(o1$coef))
		}
	}
	if(eqlabelsskip) latex_add_line(file1," \\\\")
	latex_add_line(file1,"\\textbf{Independent Variables:}\\\\")
	for(j in 1:m) {
		if(length(sigbase) == 0) {
			latex_table_entry(file1,paste("\\hspace{.4cm}",names[j],sep=""),as.coef(o1$coef[j,],o1$se[j,],k=kcoef))
		} else if(length(sigbase) == n) {
			latex_table_entry(file1,paste("\\hspace{.4cm}",names[j],sep=""),str_replace_all(as.coef(o1$coef[j,],o1$se[j,],k=kcoef,sigbase=sigbase),"\\^","\\\\^{}"))
		} else {
			stop("error in latex_model -- incompadible dimensions")
		}
		latex_table_entry(file1,"",as.se(o1$se[j,],k=kcoef))
	}
	latex_add_line(file1,"\\\\")

	if(length(stats) > 0) {
		for(j in 1:length(stats1)) {
			dig_stats <- as.matrix(o1$stats[stats1,])[j,]
			if(length(kstat) == 1) {
				if(rownames(o1$stats)[j] == "N" | rownames(o1$stats)[j] == "Clusters" | rownames(o1$stats)[j] == "C") {
					dig_stats <- dig(dig_stats,0)
				} else {
					if(length(kstat) == 1) {
						dig_stats <- dig(dig_stats,k=kstat)
					} else {
						dig_stats <- dig(dig_stats,k=kstat[j])
					}
				}
			} else {
				dig_stats <- dig(dig_stats,k=kstat[j])
			}
			if(rownames(o1$stats)[stats1[j]] == "r.squared") {
				latex_table_entry(file1,"\\textbf{R-Squared}",dig_stats)
			} else {
				latex_table_entry(file1,paste("\\textbf{",rownames(o1$stats)[stats1[j]],"}",sep=""),dig_stats)
			}
		}
	}
}

block_diag <- function(Vs)
{
	V <- Vs[[1]]
	n <- length(Vs)
	for(i in 2:n)
	{
		m1 <- nrow(Vs[[i]])
		m2 <- nrow(V)
		zmat <- matrix(rep(0,m1*m2),nrow=m1)
		V <- rbind(cbind(V,t(zmat)),cbind(zmat,Vs[[i]]))
	}
	V
}

# mediation function
func11 <- function(x) { c(x[2] +  x[5] * x[9],x[3] + x[6] * x[11],x[4] + x[7] * x[13]) }
func21 <- function(x) { c(x[2],x[3],x[4]) }
func31 <- function(x) { c(x[5] * x[9],x[6] * x[11],x[7] * x[13]) }
func41 <- function(x) { c((x[5] * x[9]) / (x[2] +  x[5] * x[9]),(x[6] * x[11]) / (x[3] + x[6] * x[11]),(x[7] * x[13]) / (x[4] + x[7] * x[13])) }

mediation1 <- function(m1,m2,m3,m4)
{
	coef <- c(m1$coef[1:7],m2$coef[1:2],m3$coef[1:2],m4$coef[1:2])
	V <- block_diag(list(V1=m1$V[1:7,1:7],V2g=m2$V[1:2,1:2],V2u=m3$V[1:2,1:2],V2i=m4$V[1:2,1:2]))

	dmt <- delta.method(func11,coef,V)
	dmd <- delta.method(func21,coef,V)
	dmi <- delta.method(func31,coef,V)
	dmr <- delta.method(func41,coef,V)

	res <- matrix(rep(NA,8*3),nrow=8)
	for(i in 1:3)
	{
		res[1,i] <- as.coef(dmt$est[i],dmt$se[i])
		res[2,i] <- as.se(dmt$se[i])
		res[3,i] <- as.coef(dmd$est[i],dmd$se[i])
		res[4,i] <- as.se(dmd$se[i])
		res[5,i] <- as.coef(dmi$est[i],dmi$se[i])
		res[6,i] <- as.se(dmi$se[i])
		res[7,i] <- dig(dmr$est[i])
		res[8,i] <- as.se(dmr$se[i])
	}

	colnames(res) <- c("Growth","Unemployment","Inflation")
	rownames(res) <- c("Total Effect","","Direct Effect","","Indirect Effect","","Indirect Ratio","")

	res
}

# mediation function, no inflation
func112 <- function(x) { c(x[2] +  x[4] * x[7],x[3] + x[5] * x[9]) }
func212 <- function(x) { c(x[2],x[3]) }
func312 <- function(x) { c(x[4] * x[7],x[5] * x[9]) }
func412 <- function(x) { c((x[4] * x[7]) / (x[2] +  x[4] * x[7]),(x[5] * x[9]) / (x[3] + x[5] * x[9])) }

mediation2 <- function(m1,m2,m3)
{
	coef <- c(m1$coef[1:5],m2$coef[1:2],m3$coef[1:2])
	V <- block_diag(list(V1=m1$V[1:5,1:5],V2g=m2$V[1:2,1:2],V2u=m3$V[1:2,1:2]))

	dmt <- delta.method(func112,coef,V)
	dmd <- delta.method(func212,coef,V)
	dmi <- delta.method(func312,coef,V)
	dmr <- delta.method(func412,coef,V)

	res <- matrix(rep(NA,8*2),nrow=8)
	for(i in 1:2)
	{
		res[1,i] <- as.coef(dmt$est[i],dmt$se[i])
		res[2,i] <- as.se(dmt$se[i])
		res[3,i] <- as.coef(dmd$est[i],dmd$se[i])
		res[4,i] <- as.se(dmd$se[i])
		res[5,i] <- as.coef(dmi$est[i],dmi$se[i])
		res[6,i] <- as.se(dmi$se[i])
		res[7,i] <- dig(dmr$est[i])
		res[8,i] <- as.se(dmr$se[i])
	}

	colnames(res) <- c("Growth","Unemployment")
	rownames(res) <- c("Total Effect","","Direct Effect","","Indirect Effect","","Indirect Ratio","")

	res
}

effsizes <- function(x,sd)
{
	gr_sd <- sd[1]
	un_sd <- sd[2]
	inf_sd <- sd[3]
	res <- rep(NA,20)

	# direct effects
	res[1] <- x[2] * gr_sd
	res[2] <- -x[3] * un_sd
	res[3] <- -x[4] * inf_sd
	res[4] <- res[1] + res[2] + res[3]

	# indirect effects
	res[5] <- x[9]  * x[5] * gr_sd
	res[6] <- x[11] * x[6] * un_sd
	res[7] <- x[13] * x[7] * inf_sd
	res[8] <- res[5] + res[6] + res[7]

	# media effects
	res[9] <- x[5] * sd[4]
	res[10] <- x[6] * sd[5]
	res[11] <- x[7] * sd[6]
	res[12] <- res[9] + res[10] + res[11]

	# total	
	res[13] <- res[1] + res[5]
	res[14] <- res[2] + res[6]
	res[15] <- res[3] + res[7]
	res[16] <- res[13] + res[14] + res[15]
	
	# direct + media
	res[17] <- res[1] + res[9]
	res[18] <- res[2] + res[10]
	res[19] <- res[3] + res[11]
	res[20] <- res[17] + res[18] + res[19]
	res
}

effsizes2 <- function(x,ch) { c(x[5] * ch[1] + x[6] * ch[2] + x[7] * ch[3]) }

#################
# load the data #
#################

data1 <- read.csv(paste(dir1,"econ_mediation_data.csv",sep=""))
country <- data1$country
monthn <- data1$monthn
year <- data1$year
newsid <- data1$newsid # unique newspaper id
monthid <- data1$monthid # unique time id

# political variables
votelead <- data1$votelead
voteleadpoll <- data1$voteleadpoll
votegov <- data1$votegov
votegovpoll <- data1$votegovpoll

# the economy
gr <- data1$gr
un <- data1$un
inf <- data1$inf

# sentiment
econ_frac_pos <- data1$econ_frac_pos
gr_frac_pos <- data1$gr_frac_pos
un_frac_pos <- data1$un_frac_pos
inf_frac_pos <- data1$inf_frac_pos

# covariates
proptype <- data1$proptype
mixedtype <- data1$mixedtype
coalsize <- data1$coalsize
timelead <- data1$timelead

# left-wing paper
first <- 1 * (newsid %% 2 == 1)

# vote intent
voteleadint <- votelead
voteleadint[is.na(voteleadint)] <- voteleadpoll[is.na(voteleadint)]
votegovint <- votegov
votegovint[is.na(votegovint)] <- votegovpoll[is.na(votegovint)]

# lagged political variables
m <- nrow(data1)
votelead_m1 <- rep(NA,m)
for(i in 2:m) {
	if(!is.na(votelead[i])) {
		for(j in (i-1):1) {
			if(newsid[j] != newsid[i]) break
			if(!is.na(votelead[j])) {
				votelead_m1[i] <- votelead[j]
				break
			}
		}
	}
}

voteleadint_m1 <- rep(NA,m)
for(i in 2:m) {
	if(!is.na(voteleadint[i])) {
		for(j in (i-1):1) {
			if(newsid[j] != newsid[i]) break
			if(!is.na(voteleadint[j])) {
				voteleadint_m1[i] <- voteleadint[j]
				break
			}
		}
	}
}

# yearly economic data
gr_yr <- 100*(exp(rowSums(log(1+panel.lag(gr,monthid,newsid,lags=0:11)/100)))-1)
un_yr <- rowMeans(panel.lag(un,monthid,newsid,lags=0:11))
unch_yr <- as.vector(un_yr - panel.lag(un_yr,monthid,newsid,lags=12))
inf_yr <- 100*(exp(rowSums(log(1+panel.lag(inf,monthid,newsid,lags=0:11)/100)))-1)
infch_yr <- as.vector(inf_yr - panel.lag(inf_yr,monthid,newsid,lags=12))

econ_frac_pos_oth <- econ_frac_pos[fmatch(paste(country,monthid,1-first,sep="_"),paste(country,monthid,first,sep="_"))]
gr_frac_pos_oth <- gr_frac_pos[fmatch(paste(country,monthid,1-first,sep="_"),paste(country,monthid,first,sep="_"))]
un_frac_pos_oth <- un_frac_pos[fmatch(paste(country,monthid,1-first,sep="_"),paste(country,monthid,first,sep="_"))]
inf_frac_pos_oth <- inf_frac_pos[fmatch(paste(country,monthid,1-first,sep="_"),paste(country,monthid,first,sep="_"))]

gr_frac_pos_avg <- rowMeans(cbind(gr_frac_pos,gr_frac_pos_oth),na.rm=T)
un_frac_pos_avg <- rowMeans(cbind(un_frac_pos,un_frac_pos_oth),na.rm=T)
inf_frac_pos_avg <- rowMeans(cbind(inf_frac_pos,inf_frac_pos_oth),na.rm=T)

###########################
# FIGURE 1 - Data Summary #
###########################

postscript(file=paste(dir1,"newsts.eps",sep=""),width=14,height=7)
countries2 <- c("Australia","Austria","Canada","France","Germany","Ireland","Israel","Italy","Japan","Luxembourg","New Zealand","Portugal","Spain","Switzerland","United Kingdom","United States")
monthyear <- paste(data1$monthn,substr(lapply(year,toString),3,4),sep="/")
par(oma=c(0,0,3,0))
par(mar=c(2,2,2,2))
par(mfrow=c(4,4))
for(i in 1:16)
{
	if(i != 14) {
		leftnews <- 2 * i - 1
		rightnews <- 2 * i
		minfrac <- min(subset(econ_frac_pos,newsid==leftnews|newsid==rightnews),na.rm=TRUE)
		maxfrac <- max(subset(econ_frac_pos,newsid==leftnews|newsid==rightnews),na.rm=TRUE)
		minmonth <- min(min(subset(monthid,newsid==leftnews&!is.na(econ_frac_pos))),min(subset(monthid,newsid==rightnews&!is.na(econ_frac_pos))))
		maxmonth <- max(max(subset(monthid,newsid==leftnews&!is.na(econ_frac_pos))),max(subset(monthid,newsid==rightnews&!is.na(econ_frac_pos))))
		diffmonth <- ceiling(.2 * (maxmonth - minmonth))
		at1 <- c(minmonth,minmonth+diffmonth,minmonth+2*diffmonth,minmonth+3*diffmonth,minmonth+4*diffmonth,minmonth+5*diffmonth)
		labels1 <- monthyear[at1]
		plot(0:1,type='n',ylim=c(minfrac-0.03,maxfrac+0.03),xlim=c(minmonth,maxmonth),ylab='Econ. Sent.',xlab='Date',main=countries2[i],at=at1,labels=labels1,cex.axis=0.8)
		points(subset(monthid,newsid==leftnews),subset(econ_frac_pos,newsid==leftnews),type='l',col='red')
		points(subset(monthid,newsid==rightnews),subset(econ_frac_pos,newsid==rightnews),type='l',col='blue',lty=2)
	}
}
dev.off()

##############################
# FIGURE 2 - Check Poll Data #
##############################

postscript(file=paste(dir1,"pollfig.eps",sep=""),width=6,height=6)
plot(votelead,voteleadpoll,xlab="Election Outcome",ylab="Poll Estimate")
abline(0, 1)
dev.off()

#########################
# Table 1 - Main Models # 
#########################

# tables 1, 2, and 3 (dropping switzerland leads to some changes)
lm1a <- lm0(votelead ~ gr_yr + un_yr + inf_yr,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1b <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1c <- lm0(votelead ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1d <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1e <- ivreg0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos + un_frac_pos + inf_frac_pos | gr_yr + un_yr + inf_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth,se="robust",timeid=monthid,unitid=country)

wald.test(function(x) x[5:6],lm1d$coef,lm1d$V)
lm1d$coef[5]
lm1d$coef[6]

tab1 <- matrix(rep("",5*17),nrow=17)

tab1[1,1] <- as.coef(lm1a$coef[1],lm1a$se[1])
tab1[2,1] <- as.se(lm1a$se[1])
tab1[3,1] <- as.coef(lm1a$coef[2],lm1a$se[2])
tab1[4,1] <- as.se(lm1a$se[2])
tab1[5,1] <- as.coef(lm1a$coef[3],lm1a$se[3])
tab1[6,1] <- as.se(lm1a$se[3])
tab1[7,1] <- as.coef(lm1a$coef[4],lm1a$se[4])
tab1[8,1] <- as.se(lm1a$se[4])
tab1[15,1] <- lm1a$N
tab1[16,1] <- lm1a$Clusters
tab1[17,1] <- dig3(lm1a$r.squared)

tab1[1,2] <- as.coef(lm1b$coef[1],lm1b$se[1])
tab1[2,2] <- as.se(lm1b$se[1])
tab1[3,2] <- as.coef(lm1b$coef[2],lm1b$se[2])
tab1[4,2] <- as.se(lm1b$se[2])
tab1[5,2] <- as.coef(lm1b$coef[3],lm1b$se[3])
tab1[6,2] <- as.se(lm1b$se[3])
tab1[7,2] <- as.coef(lm1b$coef[4],lm1b$se[4])
tab1[8,2] <- as.se(lm1b$se[4])
tab1[15,2] <- lm1b$N
tab1[16,2] <- lm1b$Clusters
tab1[17,2] <- dig3(lm1b$r.squared)

tab1[1,3] <- as.coef(lm1c$coef[1],lm1c$se[1])
tab1[2,3] <- as.se(lm1c$se[1])
tab1[3,3] <- as.coef(lm1c$coef[2],lm1c$se[2])
tab1[4,3] <- as.se(lm1c$se[2])
tab1[5,3] <- as.coef(lm1c$coef[3],lm1c$se[3])
tab1[6,3] <- as.se(lm1c$se[3])
tab1[7,3] <- as.coef(lm1c$coef[4],lm1c$se[4])
tab1[8,3] <- as.se(lm1c$se[4])
tab1[9,3] <- as.coef(lm1c$coef[5],lm1c$se[5])
tab1[10,3] <- as.se(lm1c$se[5])
tab1[11,3] <- as.coef(lm1c$coef[6],lm1c$se[6])
tab1[12,3] <- as.se(lm1c$se[6])
tab1[13,3] <- as.coef(lm1c$coef[7],lm1c$se[7])
tab1[14,3] <- as.se(lm1c$se[7])
tab1[15,3] <- lm1c$N
tab1[16,3] <- lm1c$Clusters
tab1[17,3] <- dig3(lm1c$r.squared)

tab1[1,4] <- as.coef(lm1d$coef[1],lm1d$se[1])
tab1[2,4] <- as.se(lm1d$se[1])
tab1[3,4] <- as.coef(lm1d$coef[2],lm1d$se[2])
tab1[4,4] <- as.se(lm1d$se[2])
tab1[5,4] <- as.coef(lm1d$coef[3],lm1d$se[3])
tab1[6,4] <- as.se(lm1d$se[3])
tab1[7,4] <- as.coef(lm1d$coef[4],lm1d$se[4])
tab1[8,4] <- as.se(lm1d$se[4])
tab1[9,4] <- as.coef(lm1d$coef[5],lm1d$se[5])
tab1[10,4] <- as.se(lm1d$se[5])
tab1[11,4] <- as.coef(lm1d$coef[6],lm1d$se[6])
tab1[12,4] <- as.se(lm1d$se[6])
tab1[13,4] <- as.coef(lm1d$coef[7],lm1d$se[7])
tab1[14,4] <- as.se(lm1d$se[7])
tab1[15,4] <- lm1d$N
tab1[16,4] <- lm1d$Clusters
tab1[17,4] <- dig3(lm1d$r.squared)

# econ and media, correct for measurement error
tab1[1,5] <- as.coef(lm1e$coef[1],lm1e$se[1])
tab1[2,5] <- as.se(lm1e$se[1])
tab1[3,5] <- as.coef(lm1e$coef[2],lm1e$se[2])
tab1[4,5] <- as.se(lm1e$se[2])
tab1[5,5] <- as.coef(lm1e$coef[3],lm1e$se[3])
tab1[6,5] <- as.se(lm1e$se[3])
tab1[7,5] <- as.coef(lm1e$coef[4],lm1e$se[4])
tab1[8,5] <- as.se(lm1e$se[4])
tab1[9,5] <- as.coef(lm1e$coef[5],lm1e$se[5])
tab1[10,5] <- as.se(lm1e$se[5])
tab1[11,5] <- as.coef(lm1e$coef[6],lm1e$se[6])
tab1[12,5] <- as.se(lm1e$se[6])
tab1[13,5] <- as.coef(lm1e$coef[7],lm1e$se[7])
tab1[14,5] <- as.se(lm1e$se[7])
tab1[15,5] <- lm1e$N
tab1[16,5] <- lm1e$Clusters
tab1[17,5] <- NA

tab1

#######################
# Table 2 - Mediation #
#######################

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,se="robust",unitid=country,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,se="robust",unitid=country,timeid=monthid)
lm2uc <- lm0(un_frac_pos  ~ unch_yr,se="robust",unitid=country,timeid=monthid)
lm2i  <- lm0(inf_frac_pos ~ inf_yr,se="robust",unitid=country,timeid=monthid)
lm2ic <- lm0(inf_frac_pos ~ infch_yr,se="robust",unitid=country,timeid=monthid)

med1 <- mediation1(lm1c,lm2g,lm2u,lm2i)
med2 <- mediation1(lm1d,lm2g,lm2u,lm2i)
med3 <- mediation1(lm1e,lm2g,lm2u,lm2i)

tab2 <- matrix(rep(NA,24*3),nrow=24)

tab2[1:8,1] <- med1[1:8,1]
tab2[9:16,1] <- med1[1:8,2]
tab2[17:24,1] <- med1[1:8,3]

tab2[1:8,2] <- med2[1:8,1]
tab2[9:16,2] <- med2[1:8,2]
tab2[17:24,2] <- med2[1:8,3]

tab2[1:8,3] <- med3[1:8,1]
tab2[9:16,3] <- med3[1:8,2]
tab2[17:24,3] <- med3[1:8,3]

tab2

###########################
# Table 3 -- effect sizes #
###########################

lm1e <- ivreg0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos + un_frac_pos + inf_frac_pos | gr_yr + un_yr + inf_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth,se="robust",timeid=monthid,unitid=country)

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,se="robust",unitid=country,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,se="robust",unitid=country,timeid=monthid)
lm2i  <- lm0(inf_frac_pos ~ inf_yr,se="robust",unitid=country,timeid=monthid)

sd <- c(sd0(gr_yr),sd0(un_yr),sd0(remove.outliers(inf_yr,cutoff=5)),sd0(gr_frac_pos_avg),sd0(un_frac_pos_avg),sd0(inf_frac_pos_avg))
d <- delta.method(effsizes,c(lm1e$coef,lm2g$coef,lm2u$coef,lm2i$coef),block_diag(list(V1=lm1e$V,V2g=lm2g$V,V2u=lm2u$V,V2i=lm2i$V)),sd=sd)

tab3 <- matrix(rep(NA,4*10),ncol=4)
tab3[1,] <- t(matrix(as.coef(d$est,d$se),ncol=5))[1,]
tab3[2,] <- t(matrix(as.se(d$se),ncol=5))[1,]
tab3[3,] <- t(matrix(as.coef(d$est,d$se),ncol=5))[2,]
tab3[4,] <- t(matrix(as.se(d$se),ncol=5))[2,]
tab3[5,] <- t(matrix(as.coef(d$est,d$se),ncol=5))[3,]
tab3[6,] <- t(matrix(as.se(d$se),ncol=5))[3,]
tab3[7,] <- t(matrix(as.coef(d$est,d$se),ncol=5))[4,]
tab3[8,] <- t(matrix(as.se(d$se),ncol=5))[4,]
tab3[9,] <- t(matrix(as.coef(d$est,d$se),ncol=5))[5,]
tab3[10,] <- t(matrix(as.se(d$se),ncol=5))[5,]
tab3

#############################
# Figure 3 - Aus. Sentiment #
#############################

date1 <- as.Date(paste(data1$monthn,"1",year,sep="/"),"%m/%d/%Y")
datenum <- function(x) { as.numeric(as.Date(x,"%m/%d/%Y")) }
gr_age_norm <- (gr_frac_pos - mean0(gr_frac_pos,subset=newsid==1)) / sd0(gr_frac_pos,subset=newsid==1)
gr_her_norm <- (gr_frac_pos - mean0(gr_frac_pos,subset=newsid==2)) / sd0(gr_frac_pos,subset=newsid==2)
gr_norm <- (gr - mean0(gr,subset=newsid==31)) / sd0(gr,subset=newsid==31)

postscript(file=paste(dir1,"aus.eps",sep=""),width=8,height=4)
plot(0:1,type='n',ylim=c(-3,3),xlim=c(8340-30*24+365*18-90,8340+30*24+365*18-90),xlab="",ylab="Standardized Growth and Growth Sentiment",xaxt='n',cex=.5,cex.lab=.75)
points0(date1,gr_her_norm,subset=newsid==2,type='l',col='blue',lty=3)
points0(date1,gr_age_norm,subset=newsid==1,type='l',col='red',lty=2)
points0(date1,gr_norm,subset=newsid==1,type='l')
abline(v=14842,lty=1)
adates <- c("8/1/2008","8/1/2009","8/1/2010","8/1/2011","8/1/2012")
axis(1,at=datenum(adates),labels=adates)
legend(x=c(7600+365*18-90,14650),y=c(1.25,3),c("Growth","The Herald Sun","The Age"),col=c('black','red','blue'),lty=c(1,2,3),cex=.75)
text(7600+365*18-90+860,-2.9,"Election Date",cex=.75)
dev.off()

#############################
# Figure 4 - Bush Sentiment #
#############################

date1 <- as.Date(paste(data1$monthn,"1",year,sep="/"),"%m/%d/%Y")
datenum <- function(x) { as.numeric(as.Date(x,"%m/%d/%Y")) }
gr_nyt_norm <- (gr_frac_pos - mean0(gr_frac_pos,subset=newsid==31)) / sd0(gr_frac_pos,subset=newsid==31)
gr_wsj_norm <- (gr_frac_pos - mean0(gr_frac_pos,subset=newsid==32)) / sd0(gr_frac_pos,subset=newsid==32)
gr_norm <- (gr - mean0(gr,subset=newsid==31)) / sd0(gr,subset=newsid==31)

postscript(file=paste(dir1,"bush.eps",sep=""),width=8,height=4)
plot(0:1,type='n',ylim=c(-3,3),xlim=c(8340-30*24,8340+30*24),xlab="",ylab="Standardized Growth and Growth Sentiment",xaxt='n',cex.lab=.75)
points0(date1,gr_nyt_norm,subset=newsid==31,type='l',col='red',lty=3)
points0(date1,gr_wsj_norm,subset=newsid==32,type='l',col='blue',lty=2)
points0(date1,gr_norm,subset=newsid==31,type='l')
abline(v=8340,lty=1)
adates <- c("11/1/1990","11/1/1991","11/1/1992","11/1/1993","11/1/1994")
axis(1,at=datenum(adates),labels=adates)
text(7600+840,-2.9,"Election Date",cex=.75)
legend(x=c(7600,8200),y=c(1.2,3),c("Growth","New York Times","Wall Street Journal"),col=c('black','red','blue'),lty=c(1,2,3),cex=.75)
dev.off()

### appendix tables ###

#######################################
# Table A.2 -- sequential g-estimation #
#######################################

elecsamp1     <- newsid <= 32 & first == 1 & country != "Switzerland"  & !is.na(votelead) & !is.na(gr_frac_pos_avg) & !is.na(un_frac_pos_avg) & !is.na(inf_frac_pos_avg)
elecsamp2     <- newsid <= 32 & first == 1 & country != "Switzerland"  & !is.na(voteleadint) & !is.na(gr_frac_pos_avg) & !is.na(un_frac_pos_avg) & !is.na(inf_frac_pos_avg)
df1 <- data.frame(votelead,gr_yr,un_yr,inf_yr,gr_frac_pos_avg,un_frac_pos_avg,inf_frac_pos_avg)[elecsamp1,]
df2 <- data.frame(voteleadint,gr_yr,un_yr,inf_yr,gr_frac_pos_avg,un_frac_pos_avg,inf_frac_pos_avg)[elecsamp2,]
f1 <- votelead ~ gr_yr + un_yr + inf_yr | 1 | gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg
f2 <- voteleadint ~ gr_yr + un_yr + inf_yr | 1 | gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg
de1 <- sequential_g(f1,df1)
de2 <- sequential_g(f2,df2)

# bootstrap
set.seed(1234)
R <- 100
de1boot <- matrix(rep(NA,3*100),ncol=3)
de2boot <- matrix(rep(NA,3*100),ncol=3)
for(r in 1:100) {
	print0("r",r)
	n1 <- nrow(df1)
	df1b <- df1[sample(1:n1,replace=T),]
	n2 <- nrow(df2)
	df2b <- df2[sample(1:n2,replace=T),]
	de1boot[r,] <- sequential_g(df1,df1b)$coef[2:4]
	de2boot[r,] <- sequential_g(df2,df2b)$coef[2:4]
}

tabA2 <- matrix(rep(NA,6*2),nrow=6)

tabA2[1,1] <- as.coef(de1$coef[2],sd0(de1boot[,1]))
tabA2[1,2] <- as.coef(de2$coef[2],sd0(de2boot[,1]))
tabA2[3,1] <- as.coef(de1$coef[3],sd0(de1boot[,2]))
tabA2[3,2] <- as.coef(de2$coef[3],sd0(de2boot[,2]))
tabA2[5,1] <- as.coef(de1$coef[4],sd0(de1boot[,3]))
tabA2[5,2] <- as.coef(de2$coef[4],sd0(de2boot[,3]))

tabA2[2,1] <- as.se(sd0(de1boot[,1]))
tabA2[2,2] <- as.se(sd0(de2boot[,1]))
tabA2[4,1] <- as.se(sd0(de1boot[,2]))
tabA2[4,2] <- as.se(sd0(de2boot[,2]))
tabA2[6,1] <- as.se(sd0(de1boot[,3]))
tabA2[6,2] <- as.se(sd0(de2boot[,3]))

tabA2

##########################
# Table A.3 - Media Models #
##########################

v1a <- vuongtest.lm(un_frac_pos   ~ un_yr | unch_yr)
v1b <- vuongtest.lm(inf_frac_pos  ~ inf_yr | infch_yr)

tabA3 <- matrix(rep(NA,17*5),ncol=5)

tabA3[1:4,1]   <- alternate(as.coef(lm2g$coef,lm2g$se),as.se(lm2g$se))
tabA3[1:2,2]   <- c(as.coef(lm2u$coef[1],lm2u$se[1]),as.se(lm2u$se[1]))
tabA3[5:6,2]   <- c(as.coef(lm2u$coef[2],lm2u$se[2]),as.se(lm2u$se[2]))
tabA3[1:2,3]   <- c(as.coef(lm2uc$coef[1],lm2uc$se[1]),as.se(lm2uc$se[1]))
tabA3[7:8,3]   <- c(as.coef(lm2uc$coef[2],lm2uc$se[2]),as.se(lm2uc$se[2]))
tabA3[1:2,4]   <- c(as.coef(lm2i$coef[1],lm2i$se[1]),as.se(lm2i$se[1]))
tabA3[9:10,4]  <- c(as.coef(lm2i$coef[2],lm2i$se[2]),as.se(lm2i$se[2]))
tabA3[1:2,5]   <- c(as.coef(lm2ic$coef[1],lm2ic$se[1]),as.se(lm2ic$se[1]))
tabA3[11:12,5] <- c(as.coef(lm2ic$coef[2],lm2ic$se[2]),as.se(lm2ic$se[2]))

tabA3[13,1] <- lm2g$N
tabA3[13,2] <- lm2u$N
tabA3[13,3] <- lm2uc$N
tabA3[13,4] <- lm2i$N
tabA3[13,5] <- lm2ic$N

tabA3[14,1] <- lm2g$Clusters
tabA3[14,2] <- lm2u$Clusters
tabA3[14,3] <- lm2uc$Clusters
tabA3[14,4] <- lm2i$Clusters
tabA3[14,5] <- lm2ic$Clusters

tabA3[15,1] <- dig(lm2g$r.squared)
tabA3[15,2] <- dig(lm2u$r.squared)
tabA3[15,3] <- dig(lm2uc$r.squared)
tabA3[15,4] <- dig(lm2i$r.squared)
tabA3[15,5] <- dig(lm2ic$r.squared)

tabA3[16,2] <- pval.as.pval(v1a$pval1)
tabA3[16,4] <- pval.as.pval(v1b$pval1)

tabA3[17,2] <- pval.as.pval(v1a$pval2)
tabA3[17,4] <- pval.as.pval(v1b$pval2)

tabA3

###########################################
# Table A.4 - Models (changes in inflation) #
###########################################

lm1a <- lm0(voteleadint ~ gr_yr + un_yr + infch_yr,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1b <- lm0(voteleadint ~ gr_yr + un_yr + infch_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1c <- ivreg0(voteleadint ~ gr_yr + un_yr + infch_yr + gr_frac_pos + un_frac_pos + inf_frac_pos | gr_yr + un_yr + infch_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth,se="robust",timeid=monthid,unitid=country)

v1a <- vuongtest.lm(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg | gr_yr + un_yr + infch_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1)

tabA4 <- matrix(rep("",3*19),ncol=3)

tabA4[1,1] <- as.coef(lm1a$coef[1],lm1a$se[1])
tabA4[2,1] <- as.se(lm1a$se[1])
tabA4[3,1] <- as.coef(lm1a$coef[2],lm1a$se[2])
tabA4[4,1] <- as.se(lm1a$se[2])
tabA4[5,1] <- as.coef(lm1a$coef[3],lm1a$se[3])
tabA4[6,1] <- as.se(lm1a$se[3])
tabA4[7,1] <- as.coef(lm1a$coef[4],lm1a$se[4])
tabA4[8,1] <- as.se(lm1a$se[4])
tabA4[15,1] <- lm1a$N
tabA4[16,1] <- lm1a$Clusters
tabA4[17,1] <- dig3(lm1a$r.squared)

tabA4[1,2] <- as.coef(lm1b$coef[1],lm1b$se[1])
tabA4[2,2] <- as.se(lm1b$se[1])
tabA4[3,2] <- as.coef(lm1b$coef[2],lm1b$se[2])
tabA4[4,2] <- as.se(lm1b$se[2])
tabA4[5,2] <- as.coef(lm1b$coef[3],lm1b$se[3])
tabA4[6,2] <- as.se(lm1b$se[3])
tabA4[7,2] <- as.coef(lm1b$coef[4],lm1b$se[4])
tabA4[8,2] <- as.se(lm1b$se[4])
tabA4[9,2] <- as.coef(lm1b$coef[5],lm1b$se[5])
tabA4[10,2] <- as.se(lm1b$se[5])
tabA4[11,2] <- as.coef(lm1b$coef[6],lm1b$se[6])
tabA4[12,2] <- as.se(lm1b$se[6])
tabA4[13,2] <- as.coef(lm1b$coef[7],lm1b$se[7])
tabA4[14,2] <- as.se(lm1b$se[7])
tabA4[15,2] <- lm1b$N
tabA4[16,2] <- lm1b$Clusters
tabA4[17,2] <- dig3(lm1b$r.squared)

tabA4[18,2] <- pval.as.pval(v1a$pval1)
tabA4[19,2] <- pval.as.pval(v1a$pval2)

tabA4[1,3] <- as.coef(lm1c$coef[1],lm1c$se[1])
tabA4[2,3] <- as.se(lm1c$se[1])
tabA4[3,3] <- as.coef(lm1c$coef[2],lm1c$se[2])
tabA4[4,3] <- as.se(lm1c$se[2])
tabA4[5,3] <- as.coef(lm1c$coef[3],lm1c$se[3])
tabA4[6,3] <- as.se(lm1c$se[3])
tabA4[7,3] <- as.coef(lm1c$coef[4],lm1c$se[4])
tabA4[8,3] <- as.se(lm1c$se[4])
tabA4[9,3] <- as.coef(lm1c$coef[5],lm1c$se[5])
tabA4[10,3] <- as.se(lm1c$se[5])
tabA4[11,3] <- as.coef(lm1c$coef[6],lm1c$se[6])
tabA4[12,3] <- as.se(lm1c$se[6])
tabA4[13,3] <- as.coef(lm1c$coef[7],lm1c$se[7])
tabA4[14,3] <- as.se(lm1c$se[7])
tabA4[15,3] <- lm1c$N
tabA4[16,3] <- lm1c$Clusters
tabA4[17,3] <- NA

tabA4

#############################################
# Table A.5 - Mediation, Changes in Inflation #
#############################################

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,se="robust",unitid=country,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,se="robust",unitid=country,timeid=monthid)
lm2uc <- lm0(un_frac_pos  ~ unch_yr,se="robust",unitid=country,timeid=monthid)
lm2i  <- lm0(inf_frac_pos ~ inf_yr,se="robust",unitid=country,timeid=monthid)
lm2ic <- lm0(inf_frac_pos ~ infch_yr,se="robust",unitid=country,timeid=monthid)

med1 <- mediation1(lm1b,lm2g,lm2u,lm2ic)
med2 <- mediation1(lm1c,lm2g,lm2u,lm2ic)

tabA5 <- matrix(rep(NA,24*2),nrow=24)

tabA5[1:8,1] <- med1[1:8,1]
tabA5[9:16,1] <- med1[1:8,2]
tabA5[17:24,1] <- med1[1:8,3]

tabA5[1:8,2] <- med2[1:8,1]
tabA5[9:16,2] <- med2[1:8,2]
tabA5[17:24,2] <- med2[1:8,3]

tabA5

###################################################
# Table A.6 -- Mediation, Clustered Standard Errors #
###################################################

# clustered ses
lm1a <- lm0(votelead ~ gr_yr + un_yr + inf_yr,subset=first==1,se="cluster",timeid=monthid,unitid=country)
lm1b <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr,subset=first==1,se="cluster",timeid=monthid,unitid=country)
lm1c <- lm0(votelead ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="cluster",timeid=monthid,unitid=country)
lm1d <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="cluster",timeid=monthid,unitid=country)
lm1e <- ivreg0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos + un_frac_pos + inf_frac_pos | gr_yr + un_yr + inf_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth,se="cluster",timeid=monthid,unitid=country)

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,se="cluster",unitid=newsid,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,se="cluster",unitid=newsid,timeid=monthid)
lm2i  <- lm0(inf_frac_pos  ~ inf_yr,se="cluster",unitid=newsid,timeid=monthid)

med1 <- mediation1(lm1c,lm2g,lm2u,lm2i)
med2 <- mediation1(lm1d,lm2g,lm2u,lm2i)
med3 <- mediation1(lm1e,lm2g,lm2u,lm2i)

tabA6 <- matrix(rep(NA,24*3),nrow=24)

tabA6[1:8,1] <- med1[1:8,1]
tabA6[9:16,1] <- med1[1:8,2]
tabA6[17:24,1] <- med1[1:8,3]

tabA6[1:8,2] <- med2[1:8,1]
tabA6[9:16,2] <- med2[1:8,2]
tabA6[17:24,2] <- med2[1:8,3]

tabA6[1:8,3] <- med3[1:8,1]
tabA6[9:16,3] <- med3[1:8,2]
tabA6[17:24,3] <- med3[1:8,3]

tabA6

##########################################################
# Table A.7 -- Mediation, Panel Newey-West Standard Errors #
##########################################################

# clustered ses
lm1a <- lm0(votelead ~ gr_yr + un_yr + inf_yr,subset=first==1,se="pnw",timeid=monthid,unitid=country)
lm1b <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr,subset=first==1,se="pnw",timeid=monthid,unitid=country)
lm1c <- lm0(votelead ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="pnw",timeid=monthid,unitid=country)
lm1d <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="pnw",timeid=monthid,unitid=country)
lm1e <- ivreg0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos + un_frac_pos + inf_frac_pos | gr_yr + un_yr + inf_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth,se="pnw",timeid=monthid,unitid=country)

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,se="pnw",unitid=newsid,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,se="pnw",unitid=newsid,timeid=monthid)
lm2i  <- lm0(inf_frac_pos  ~ inf_yr,se="pnw",unitid=newsid,timeid=monthid)

med1 <- mediation1(lm1c,lm2g,lm2u,lm2i)
med2 <- mediation1(lm1d,lm2g,lm2u,lm2i)
med3 <- mediation1(lm1e,lm2g,lm2u,lm2i)

tabA7 <- matrix(rep(NA,24*3),nrow=24)

tabA7[1:8,1] <- med1[1:8,1]
tabA7[9:16,1] <- med1[1:8,2]
tabA7[17:24,1] <- med1[1:8,3]

tabA7[1:8,2] <- med2[1:8,1]
tabA7[9:16,2] <- med2[1:8,2]
tabA7[17:24,2] <- med2[1:8,3]

tabA7[1:8,3] <- med3[1:8,1]
tabA7[9:16,3] <- med3[1:8,2]
tabA7[17:24,3] <- med3[1:8,3]

tabA7

###############################################
# table A.8 -- level vs. changes on vote share #
###############################################

lm1  <- lm0(votelead  ~ gr_yr + un_yr + inf_yr,subset=first==1,se="robust",unitid=country,timeid=monthid)
lm2  <- lm0(votelead  ~ gr_yr + unch_yr + inf_yr,subset=first==1,se="robust",unitid=country,timeid=monthid)
lm3  <- lm0(votelead  ~ gr_yr + un_yr + infch_yr,subset=first==1,se="robust",unitid=country,timeid=monthid)
lm4  <- lm0(votelead  ~ gr_yr + unch_yr + infch_yr,subset=first==1,se="robust",unitid=country,timeid=monthid)
lm5  <- lm0(voteleadint  ~ gr_yr + un_yr + inf_yr,subset=first==1,se="robust",unitid=country,timeid=monthid)
lm6  <- lm0(voteleadint  ~ gr_yr + unch_yr + inf_yr,subset=first==1,se="robust",unitid=country,timeid=monthid)
lm7  <- lm0(voteleadint  ~ gr_yr + un_yr + infch_yr,subset=first==1,se="robust",unitid=country,timeid=monthid)
lm8  <- lm0(voteleadint  ~ gr_yr + unch_yr + infch_yr,subset=first==1,se="robust",unitid=country,timeid=monthid)

o1 <- output.models(list("(1)"=lm1,"(2)"=lm2,"(3)"=lm3,"(4)"=lm4,"(5)"=lm5,"(6)"=lm6,"(7)"=lm7,"(8)"=lm8))

v1a <- vuongtest.lm.multi(votelead    ~ gr_yr + un_yr + inf_yr | gr_yr + unch_yr + inf_yr | gr_yr + un_yr + infch_yr | gr_yr + unch_yr + infch_yr,subset=first==1)
v1b <- vuongtest.lm.multi(voteleadint ~ gr_yr + un_yr + inf_yr | gr_yr + unch_yr + inf_yr | gr_yr + un_yr + infch_yr | gr_yr + unch_yr + infch_yr,subset=first==1)

tabA8 <- matrix(rep(NA,8*15),ncol=8)
tabA8[c(1,3,5,7,9,11),] <- as.coef(o1$coef,o1$se)
tabA8[c(2,4,6,8,10,12),] <- as.se(o1$se)
tabA8[13:15,] <- o1$stats
rownames(tabA8) <- c(alternate(rownames(o1$coef),rep("",6)),rownames(o1$stats))
tabA8

###########################
# Table A.9 -- ldv version #
###########################

lm1a <- lm0(votelead ~ gr_yr + un_yr + inf_yr + votelead_m1,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1b <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr + voteleadint_m1,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1c <- lm0(votelead ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg + votelead_m1,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1d <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg + voteleadint_m1,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1e <- ivreg0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos + un_frac_pos + inf_frac_pos + voteleadint_m1 | gr_yr + un_yr + inf_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth + voteleadint_m1,se="robust",timeid=monthid,unitid=country)

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,se="robust",unitid=country,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,se="robust",unitid=country,timeid=monthid)
lm2uc <- lm0(un_frac_pos  ~ unch_yr,se="robust",unitid=country,timeid=monthid)
lm2i  <- lm0(inf_frac_pos ~ inf_yr,se="robust",unitid=country,timeid=monthid)
lm2ic <- lm0(inf_frac_pos ~ infch_yr,se="robust",unitid=country,timeid=monthid)

med1 <- mediation1(lm1c,lm2g,lm2u,lm2i)
med2 <- mediation1(lm1d,lm2g,lm2u,lm2i)
med3 <- mediation1(lm1e,lm2g,lm2u,lm2i)

tabA9 <- matrix(rep(NA,24*3),nrow=24)

tabA9[1:8,1] <- med1[1:8,1]
tabA9[9:16,1] <- med1[1:8,2]
tabA9[17:24,1] <- med1[1:8,3]

tabA9[1:8,2] <- med2[1:8,1]
tabA9[9:16,2] <- med2[1:8,2]
tabA9[17:24,2] <- med2[1:8,3]

tabA9[1:8,3] <- med3[1:8,1]
tabA9[9:16,3] <- med3[1:8,2]
tabA9[17:24,3] <- med3[1:8,3]

tabA9

### table A.10 ###

############################
# Table 10 -- No Inflation #
############################

lm1a <- lm0(votelead ~ gr_yr + un_yr,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1b <- lm0(voteleadint ~ gr_yr + un_yr,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1c <- lm0(votelead ~ gr_yr + un_yr + gr_frac_pos_avg + un_frac_pos_avg,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1d <- lm0(voteleadint ~ gr_yr + un_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1e <- ivreg0(voteleadint ~ gr_yr + un_yr + gr_frac_pos + un_frac_pos | gr_yr + un_yr + gr_frac_pos_oth + un_frac_pos_oth,se="robust",timeid=monthid,unitid=country)

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,se="robust",unitid=newsid,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,se="robust",unitid=newsid,timeid=monthid)

med1 <- mediation2(lm1c,lm2g,lm2u)
med2 <- mediation2(lm1d,lm2g,lm2u)
med3 <- mediation2(lm1e,lm2g,lm2u)

tabA10 <- matrix(rep(NA,16*3),nrow=16)

tabA10[1:8,1] <- med1[1:8,1]
tabA10[9:16,1] <- med1[1:8,2]

tabA10[1:8,2] <- med2[1:8,1]
tabA10[9:16,2] <- med2[1:8,2]

tabA10[1:8,3] <- med3[1:8,1]
tabA10[9:16,3] <- med3[1:8,2]

tabA10

###################################################
# Table A.11 -- Inflation trimmed #
###################################################

lm1a <- lm0(votelead ~ gr_yr + un_yr + inf_yr,subset=first==1&abs(inf_yr)<=10,se="robust",timeid=monthid,unitid=country)
lm1b <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr,subset=first==1&abs(inf_yr)<=10,se="robust",timeid=monthid,unitid=country)
lm1c <- lm0(votelead ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1&abs(inf_yr)<=10,se="robust",timeid=monthid,unitid=country)
lm1d <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1&abs(inf_yr)<=10,se="robust",timeid=monthid,unitid=country)
lm1e <- ivreg0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos + un_frac_pos + inf_frac_pos | gr_yr + un_yr + inf_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth,subset=abs(inf_yr)<=10,se="robust",timeid=monthid,unitid=country)

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,subset=abs(inf_yr)<=10,se="robust",unitid=newsid,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,subset=abs(inf_yr)<=10,se="robust",unitid=newsid,timeid=monthid)
lm2i  <- lm0(inf_frac_pos  ~ inf_yr,subset=abs(inf_yr)<=10,se="robust",unitid=newsid,timeid=monthid)

med1 <- mediation1(lm1c,lm2g,lm2u,lm2i)
med2 <- mediation1(lm1d,lm2g,lm2u,lm2i)
med3 <- mediation1(lm1e,lm2g,lm2u,lm2i)

tabA11 <- matrix(rep(NA,24*3),nrow=24)

tabA11[1:8,1] <- med1[1:8,1]
tabA11[9:16,1] <- med1[1:8,2]
tabA11[17:24,1] <- med1[1:8,3]

tabA11[1:8,2] <- med2[1:8,1]
tabA11[9:16,2] <- med2[1:8,2]
tabA11[17:24,2] <- med2[1:8,3]

tabA11[1:8,3] <- med3[1:8,1]
tabA11[9:16,3] <- med3[1:8,2]
tabA11[17:24,3] <- med3[1:8,3]

tabA11

###############################
# Table A.12 -- control version #
###############################

lm1a <- lm0(votelead ~ gr_yr + un_yr + inf_yr + proptype + mixedtype + coalsize + timelead,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1b <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr + proptype + mixedtype + coalsize + timelead,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1c <- lm0(votelead ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg + proptype + mixedtype + coalsize + timelead,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1d <- lm0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg + proptype + mixedtype + coalsize + timelead,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1e <- ivreg0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos + un_frac_pos + inf_frac_pos + proptype + mixedtype + coalsize + timelead | gr_yr + un_yr + inf_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth + proptype + mixedtype + coalsize + timelead,se="robust",timeid=monthid,unitid=country)

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,se="robust",unitid=country,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,se="robust",unitid=country,timeid=monthid)
lm2uc <- lm0(un_frac_pos  ~ unch_yr,se="robust",unitid=country,timeid=monthid)
lm2i  <- lm0(inf_frac_pos ~ inf_yr,se="robust",unitid=country,timeid=monthid)
lm2ic <- lm0(inf_frac_pos ~ infch_yr,se="robust",unitid=country,timeid=monthid)

med1 <- mediation1(lm1c,lm2g,lm2u,lm2i)
med2 <- mediation1(lm1d,lm2g,lm2u,lm2i)
med3 <- mediation1(lm1e,lm2g,lm2u,lm2i)

tabA12 <- matrix(rep(NA,24*3),nrow=24)

tabA12[1:8,1] <- med1[1:8,1]
tabA12[9:16,1] <- med1[1:8,2]
tabA12[17:24,1] <- med1[1:8,3]

tabA12[1:8,2] <- med2[1:8,1]
tabA12[9:16,2] <- med2[1:8,2]
tabA12[17:24,2] <- med2[1:8,3]

tabA12[1:8,3] <- med3[1:8,1]
tabA12[9:16,3] <- med3[1:8,2]
tabA12[17:24,3] <- med3[1:8,3]

tabA12

#################################
# Table A.13 -- Junior partners #
#################################

votejunior <- votegov - votelead
votejuniorint <- votegovint - voteleadint
lm1 <- lm0(votelead ~ gr_yr + un_yr + inf_yr,subset=first==1,unitid=country)
lm2 <- lm0(votejunior ~ gr_yr + un_yr + inf_yr,subset=first==1& abs(votelead - votegov)>0.001,unitid=country)
lm3 <- lm0(votegov ~ gr_yr + un_yr + inf_yr,subset=first==1,unitid=country)
tabA13 <- output.models(list("(1)"=lm1,"(2)"=lm2,"(3)"=lm3),rename=list("N"="Number of Elections","Clusters"="Number of Countries","r.squared"="R-Squared"))
print.output(tabA13)

##################################
# Table A.14 - Gov. party models # 
##################################

# tables 1, 2, and 3 (dropping switzerland leads to some changes)
lm1a <- lm0(votegov ~ gr_yr + un_yr + inf_yr,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1b <- lm0(votegovint ~ gr_yr + un_yr + inf_yr,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1c <- lm0(votegov ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1d <- lm0(votegovint ~ gr_yr + un_yr + inf_yr + gr_frac_pos_avg + un_frac_pos_avg + inf_frac_pos_avg,subset=first==1,se="robust",timeid=monthid,unitid=country)
lm1e <- ivreg0(votegovint ~ gr_yr + un_yr + inf_yr + gr_frac_pos + un_frac_pos + inf_frac_pos | gr_yr + un_yr + inf_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth,se="robust",timeid=monthid,unitid=country)

tabA14 <- matrix(rep("",5*17),nrow=17)

tabA14[1,1] <- as.coef(lm1a$coef[1],lm1a$se[1])
tabA14[2,1] <- as.se(lm1a$se[1])
tabA14[3,1] <- as.coef(lm1a$coef[2],lm1a$se[2])
tabA14[4,1] <- as.se(lm1a$se[2])
tabA14[5,1] <- as.coef(lm1a$coef[3],lm1a$se[3])
tabA14[6,1] <- as.se(lm1a$se[3])
tabA14[7,1] <- as.coef(lm1a$coef[4],lm1a$se[4])
tabA14[8,1] <- as.se(lm1a$se[4])
tabA14[15,1] <- lm1a$N
tabA14[16,1] <- lm1a$Clusters
tabA14[17,1] <- dig3(lm1a$r.squared)

tabA14[1,2] <- as.coef(lm1b$coef[1],lm1b$se[1])
tabA14[2,2] <- as.se(lm1b$se[1])
tabA14[3,2] <- as.coef(lm1b$coef[2],lm1b$se[2])
tabA14[4,2] <- as.se(lm1b$se[2])
tabA14[5,2] <- as.coef(lm1b$coef[3],lm1b$se[3])
tabA14[6,2] <- as.se(lm1b$se[3])
tabA14[7,2] <- as.coef(lm1b$coef[4],lm1b$se[4])
tabA14[8,2] <- as.se(lm1b$se[4])
tabA14[15,2] <- lm1b$N
tabA14[16,2] <- lm1b$Clusters
tabA14[17,2] <- dig3(lm1b$r.squared)

tabA14[1,3] <- as.coef(lm1c$coef[1],lm1c$se[1])
tabA14[2,3] <- as.se(lm1c$se[1])
tabA14[3,3] <- as.coef(lm1c$coef[2],lm1c$se[2])
tabA14[4,3] <- as.se(lm1c$se[2])
tabA14[5,3] <- as.coef(lm1c$coef[3],lm1c$se[3])
tabA14[6,3] <- as.se(lm1c$se[3])
tabA14[7,3] <- as.coef(lm1c$coef[4],lm1c$se[4])
tabA14[8,3] <- as.se(lm1c$se[4])
tabA14[9,3] <- as.coef(lm1c$coef[5],lm1c$se[5])
tabA14[10,3] <- as.se(lm1c$se[5])
tabA14[11,3] <- as.coef(lm1c$coef[6],lm1c$se[6])
tabA14[12,3] <- as.se(lm1c$se[6])
tabA14[13,3] <- as.coef(lm1c$coef[7],lm1c$se[7])
tabA14[14,3] <- as.se(lm1c$se[7])
tabA14[15,3] <- lm1c$N
tabA14[16,3] <- lm1c$Clusters
tabA14[17,3] <- dig3(lm1c$r.squared)

tabA14[1,4] <- as.coef(lm1d$coef[1],lm1d$se[1])
tabA14[2,4] <- as.se(lm1d$se[1])
tabA14[3,4] <- as.coef(lm1d$coef[2],lm1d$se[2])
tabA14[4,4] <- as.se(lm1d$se[2])
tabA14[5,4] <- as.coef(lm1d$coef[3],lm1d$se[3])
tabA14[6,4] <- as.se(lm1d$se[3])
tabA14[7,4] <- as.coef(lm1d$coef[4],lm1d$se[4])
tabA14[8,4] <- as.se(lm1d$se[4])
tabA14[9,4] <- as.coef(lm1d$coef[5],lm1d$se[5])
tabA14[10,4] <- as.se(lm1d$se[5])
tabA14[11,4] <- as.coef(lm1d$coef[6],lm1d$se[6])
tabA14[12,4] <- as.se(lm1d$se[6])
tabA14[13,4] <- as.coef(lm1d$coef[7],lm1d$se[7])
tabA14[14,4] <- as.se(lm1d$se[7])
tabA14[15,4] <- lm1d$N
tabA14[16,4] <- lm1d$Clusters
tabA14[17,4] <- dig3(lm1d$r.squared)

# econ and media, correct for measurement error
tabA14[1,5] <- as.coef(lm1e$coef[1],lm1e$se[1])
tabA14[2,5] <- as.se(lm1e$se[1])
tabA14[3,5] <- as.coef(lm1e$coef[2],lm1e$se[2])
tabA14[4,5] <- as.se(lm1e$se[2])
tabA14[5,5] <- as.coef(lm1e$coef[3],lm1e$se[3])
tabA14[6,5] <- as.se(lm1e$se[3])
tabA14[7,5] <- as.coef(lm1e$coef[4],lm1e$se[4])
tabA14[8,5] <- as.se(lm1e$se[4])
tabA14[9,5] <- as.coef(lm1e$coef[5],lm1e$se[5])
tabA14[10,5] <- as.se(lm1e$se[5])
tabA14[11,5] <- as.coef(lm1e$coef[6],lm1e$se[6])
tabA14[12,5] <- as.se(lm1e$se[6])
tabA14[13,5] <- as.coef(lm1e$coef[7],lm1e$se[7])
tabA14[14,5] <- as.se(lm1e$se[7])
tabA14[15,5] <- lm1e$N
tabA14[16,5] <- lm1e$Clusters
tabA14[17,5] <- NA

tabA14

#############################################
# Table A.15 - Mediation, gov. party models #
#############################################

lm2g  <- lm0(gr_frac_pos  ~ gr_yr,se="robust",unitid=country,timeid=monthid)
lm2u  <- lm0(un_frac_pos  ~ un_yr,se="robust",unitid=country,timeid=monthid)
lm2uc <- lm0(un_frac_pos  ~ unch_yr,se="robust",unitid=country,timeid=monthid)
lm2i  <- lm0(inf_frac_pos ~ inf_yr,se="robust",unitid=country,timeid=monthid)
lm2ic <- lm0(inf_frac_pos ~ infch_yr,se="robust",unitid=country,timeid=monthid)

med1 <- mediation1(lm1c,lm2g,lm2u,lm2i)
med2 <- mediation1(lm1d,lm2g,lm2u,lm2i)
med3 <- mediation1(lm1e,lm2g,lm2u,lm2i)

tabA15 <- matrix(rep(NA,24*3),nrow=24)

tabA15[1:8,1] <- med1[1:8,1]
tabA15[9:16,1] <- med1[1:8,2]
tabA15[17:24,1] <- med1[1:8,3]

tabA15[1:8,2] <- med2[1:8,1]
tabA15[9:16,2] <- med2[1:8,2]
tabA15[17:24,2] <- med2[1:8,3]

tabA15[1:8,3] <- med3[1:8,1]
tabA15[9:16,3] <- med3[1:8,2]
tabA15[17:24,3] <- med3[1:8,3]

tabA15

#############################################
# Misc. other numbers reported in the paper #
#############################################

# quality of vote intention reported in election results and vote intention section:
lm0(voteleadpoll ~ votelead)$r.squared
lm0(voteleadpoll ~ votelead)$sigma

# reliabilities reported in measurement error section:
2 * cov0(econ_frac_pos,econ_frac_pos_oth) / (sd0(econ_frac_pos)^2 + cov0(econ_frac_pos,econ_frac_pos_oth))
2 * cov0(gr_frac_pos,gr_frac_pos_oth) / (sd0(gr_frac_pos)^2 + cov0(gr_frac_pos,gr_frac_pos_oth))
2 * cov0(un_frac_pos,un_frac_pos_oth) / (sd0(un_frac_pos)^2 + cov0(un_frac_pos,un_frac_pos_oth))
2 * cov0(inf_frac_pos,inf_frac_pos_oth) / (sd0(inf_frac_pos)^2 + cov0(inf_frac_pos,inf_frac_pos_oth))

# effect sizes reported in substantive effects section:

# return the preferred model
lm1 <- ivreg0(voteleadint ~ gr_yr + un_yr + inf_yr + gr_frac_pos + un_frac_pos + inf_frac_pos | gr_yr + un_yr + inf_yr + gr_frac_pos_oth + un_frac_pos_oth + inf_frac_pos_oth,se="robust",timeid=monthid,unitid=country)

# stadardized gr
st_gr <- (gr_yr - mean0(gr_yr,subset=first==1)) / sd0(gr_yr,subset=first==1)
st_un <- (un_yr - mean0(un_yr,subset=first==1)) / sd0(un_yr,subset=first==1)
st_inf <- (inf_yr - mean0(inf_yr,subset=first==1)) / sd0(inf_yr,subset=first==1)

# stadardized sentiment
st_gr_sent <- (gr_frac_pos_avg - mean0(gr_frac_pos_avg,subset=first==1)) / sd0(gr_frac_pos_avg,subset=first==1)
st_un_sent <- (un_frac_pos_avg - mean0(un_frac_pos_avg,subset=first==1)) / sd0(un_frac_pos_avg,subset=first==1)
st_inf_sent <- (inf_frac_pos_avg - mean0(inf_frac_pos_avg,subset=first==1)) / sd0(inf_frac_pos_avg,subset=first==1)

# set standardized sentiment to standardized econ
try_gr_sent <- sd0(gr_frac_pos,subset=first==1) * st_gr + mean0(gr_frac_pos,subset=first==1)
try_un_sent <- sd0(un_frac_pos,subset=first==1) * st_un + mean0(un_frac_pos,subset=first==1)
try_inf_sent <- sd0(inf_frac_pos,subset=first==1) * st_inf + mean0(inf_frac_pos,subset=first==1)

pred1 <- cbind(1,gr_yr,un_yr,inf_yr,try_gr_sent,try_un_sent,try_inf_sent) %*% lm1$coef - cbind(1,gr_yr,un_yr,inf_yr,gr_frac_pos_avg,un_frac_pos_avg,inf_frac_pos_avg) %*% lm1$coef

mean0(pred1,subset=first==1&!is.na(voteleadint))
quantile0(pred1,subset=first==1&!is.na(voteleadint),probs=c(0.1,0.9))

# how much better media coverage was that reality (in stad. dev.), reported in saved by the media section:
samp1 <- country=="Australia"&year==2010&first==1&monthn>=2&monthn<=7
mean0(st_gr_sent[samp1] - st_gr[samp1])

# same number, in levels
(sd0(gr_frac_pos_avg,subset=first==1) * mean0(st_gr_sent[samp1] - st_gr[samp1]) + mean0(gr_frac_pos_avg,subset=first==1))

# corresponding effect on vote share
lm1$coef[5] * sd0(gr_frac_pos,subset=first==1) * (sd0(gr_frac_pos_avg,subset=first==1) * mean0(st_gr_sent[samp1] - st_gr[samp1]) + mean0(gr_frac_pos_avg,subset=first==1))

# reported in false positive section:
lm1$coef[5] * sd0(gr_frac_pos,subset=first==1)

# how much better media coverage was that reality (in stad. dev.), reported in saved by the media section:
samp1 <- country=="United States"&year==1992&first==1&monthn>=5&monthn<=10
mean0(st_gr_sent[samp1] - st_gr[samp1])

# same number, in levels
(sd0(gr_frac_pos_avg,subset=first==1) * mean0(st_gr_sent[samp1] - st_gr[samp1]) + mean0(gr_frac_pos_avg,subset=first==1))

# corresponding effect on vote share
lm1$coef[5] * sd0(gr_frac_pos,subset=first==1) * (sd0(gr_frac_pos_avg,subset=first==1) * mean0(st_gr_sent[samp1] - st_gr[samp1]) + mean0(gr_frac_pos_avg,subset=first==1))

################
# Latex Tables #
################

file1 <- paste(dir1,"econ_mediation_tables.tex",sep="")
latex_doc_start(file1)

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) & (4) & (5) \\\\")
latex_add_line(file1,"\\textbf{Dependent Variable:} & Vote Share & Vote Intent & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Independent Variables:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Constant",tab1[1,])
latex_table_entry(file1,"",tab1[2,])
latex_table_entry(file1,"\\hspace{.4cm}Growth (yearly)",tab1[3,])
latex_table_entry(file1,"",tab1[4,])
latex_table_entry(file1,"\\hspace{.4cm}Unem. (yearly)",tab1[5,])
latex_table_entry(file1,"",tab1[6,])
latex_table_entry(file1,"\\hspace{.4cm}Inf. (yearly)",tab1[7,])
latex_table_entry(file1,"",tab1[8,])
latex_table_entry(file1,"\\hspace{.4cm}Growth Sentiment",tab1[9,])
latex_table_entry(file1,"",tab1[10,])
latex_table_entry(file1,"\\hspace{.4cm}Unem. Sentiment",tab1[11,])
latex_table_entry(file1,"",tab1[12,])
latex_table_entry(file1,"\\hspace{.4cm}Inf. Sentiment",tab1[13,])
latex_table_entry(file1,"",tab1[14,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & OLS & OLS & IV \\\\")
latex_table_entry(file1,"\\textbf{Number of Months}",tab1[15,])
latex_table_entry(file1,"\\textbf{Number of Countries}",tab1[16,])
latex_table_entry(file1,"\\textbf{R-Squared}",tab1[17,])
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Elections Outcomes, the Economy, and Sentiment --- Heteroskedasticity robust standard errors in parentheses. Column (5) corrects the results of column (4) for measurement error using instrumental variables. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{VoteShCoefTable}")
latex_add_line(file1,"\\end{table}\n\n\n")

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tab2[1,])
latex_table_entry(file1,"",tab2[2,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tab2[3,])
latex_table_entry(file1,"",tab2[4,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tab2[5,])
latex_table_entry(file1,"",tab2[6,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tab2[7,])
latex_table_entry(file1,"",tab2[8,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tab2[9,])
latex_table_entry(file1,"",tab2[10,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tab2[11,])
latex_table_entry(file1,"",tab2[12,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tab2[13,])
latex_table_entry(file1,"",tab2[14,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tab2[15,])
latex_table_entry(file1,"",tab2[16,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Inflation:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tab2[17,])
latex_table_entry(file1,"",tab2[18,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tab2[19,])
latex_table_entry(file1,"",tab2[20,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tab2[21,])
latex_table_entry(file1,"",tab2[22,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tab2[23,])
latex_table_entry(file1,"",tab2[24,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & IV \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis --- Standard errors calculated using the delta method in parentheses. Results in columns (1), (2), and (3) are calculated based on columns (3), (4), and (5) of Table \\ref{VoteShCoefTable}, respectively. The first column reports mediation analysis for vote share, the second reports mediation analysis for vote intent, and the third corrects the results in column (2) for measurement error using instrumental variables. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTable}")
latex_add_line(file1,"\\end{table}\n\n\n")

# table 3
latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c c} \\hline \\hline")
latex_add_line(file1," \\textbf{Economic Variable:}  & Growth & Unem. & Inflation & All \\\\")
latex_add_line(file1,"\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tab3[1,])
latex_table_entry(file1,"",tab3[2,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tab3[3,])
latex_table_entry(file1,"",tab3[4,])
latex_table_entry(file1,"\\hspace{.4cm}Media Effect",tab3[5,])
latex_table_entry(file1,"",tab3[6,])
latex_table_entry(file1,"\\hspace{.4cm}Direct + Indirect",tab3[7,])
latex_table_entry(file1,"",tab3[8,])
latex_table_entry(file1,"\\hspace{.4cm}Direct + Media",tab3[9,])
latex_table_entry(file1,"",tab3[10,])
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Effect Sizes --- Effect sizes on intention to vote for the incumbent's party are measured for a one standard deviation improvement in growth, unemployment, inflation, and all three economic statistics simultaneously. Estimates are calculated based on column (5) of Table \\ref{VoteShCoefTable} and columns (1), (2), and (4) of Table \\ref{MediaTable}. Standard errors calculated using the delta method in parentheses. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{EffectSizeTable}")
latex_add_line(file1,"\\end{table}\n\n\n")

## TABLE A.2
latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Share & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Average Controlled Direct Effect",tabA2[1,])
latex_table_entry(file1,"",tabA2[2,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Average Controlled Direct Effect",tabA2[3,])
latex_table_entry(file1,"",tabA2[4,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Inflation:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Average Controlled Direct Effect",tabA2[5,])
latex_table_entry(file1,"",tabA2[6,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & Sequential g-Estimation & Sequential g-Estimation \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis --- Standard errors calculated using the delta method in parentheses. Results in columns (1), (2), and (3) are calculated based on columns (3), (4), and (5) of Table \\ref{VoteShCoefTable}, respectively. The first column reports mediation analysis for vote share, the second reports mediation analysis for vote intent, and the third corrects the results in column (2) for measurement error using instrumental variables. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTableGEstBoot}")
latex_add_line(file1,"\\end{table}\n\n\n")

# Table A.3

latex_add_line(file1,"\\begin{landscape}")
latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c c c } \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) & (4) & (5) \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Growth & Unemployment & Unemployment & Inflation & Inflation \\\\")
latex_add_line(file1,"  & Sentiment & Sentiment & Sentiment & Sentiment & Sentiment \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Independent Variables:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Constant",tabA3[1,])
latex_table_entry(file1,"",tabA3[2,])
latex_table_entry(file1,"\\hspace{.4cm}Growth (yearly)",tabA3[3,])
latex_table_entry(file1,"",tabA3[4,])
latex_table_entry(file1,"\\hspace{.4cm}Unemployment (yearly)",tabA3[5,])
latex_table_entry(file1,"",tabA3[6,])
latex_table_entry(file1,"\\hspace{.4cm}Unemployment Changes (yearly)",tabA3[7,])
latex_table_entry(file1,"",tabA3[8,])
latex_table_entry(file1,"\\hspace{.4cm}Inflation (yearly)",tabA3[9,])
latex_table_entry(file1,"",tabA3[10,])
latex_table_entry(file1,"\\hspace{.4cm}Inflation Changes (yearly)",tabA3[11,])
latex_table_entry(file1,"",tabA3[12,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Vuong Test of Levels vs. Changes:} \\\\")
latex_add_line(file1,paste("\\hspace{.4cm}p-Value for null that parameter is in the overlapping region & & \\multicolumn{2}{c}{",tabA3[16,2],"} & \\multicolumn{2}{c}{",tabA3[16,4] ,"} \\\\",sep=""))
latex_add_line(file1,paste("\\hspace{.4cm}p-Value for null that models fit equally well & & \\multicolumn{2}{c}{",tabA3[16,2],"} & \\multicolumn{2}{c}{",tabA3[16,4] ,"} \\\\",sep=""))
latex_add_line(file1,"\\\\")
latex_table_entry(file1,"\\textbf{Number of Months}",tabA3[13,])
latex_table_entry(file1,"\\textbf{Number of Countries}",tabA3[14,])
latex_table_entry(file1,"\\textbf{R-Squared}",tabA3[15,])
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Sentiment vs. the Economy --- Heteroskedasticity robust standard errors in parentheses. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediaTable}")
latex_add_line(file1,"\\end{table}\n\n\n")
latex_add_line(file1,"\\end{landscape}")

# table A.4 #

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c } \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3)  \\\\")
latex_add_line(file1,"\\textbf{Dependent Variable:} & Vote Intent & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Independent Variables:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Constant",tabA4[1,])
latex_table_entry(file1,"",tabA4[2,])
latex_table_entry(file1,"\\hspace{.4cm}Growth (yearly)",tabA4[3,])
latex_table_entry(file1,"",tabA4[4,])
latex_table_entry(file1,"\\hspace{.4cm}Unem. (yearly)",tabA4[5,])
latex_table_entry(file1,"",tabA4[6,])
latex_table_entry(file1,"\\hspace{.4cm}Inf. Changes (yearly)",tabA4[7,])
latex_table_entry(file1,"",tabA4[8,])
latex_table_entry(file1,"\\hspace{.4cm}Growth Sentiment",tabA4[9,])
latex_table_entry(file1,"",tabA4[10,])
latex_table_entry(file1,"\\hspace{.4cm}Unem. Sentiment",tabA4[11,])
latex_table_entry(file1,"",tabA4[12,])
latex_table_entry(file1,"\\hspace{.4cm}Inf. Sentiment",tabA4[13,])
latex_table_entry(file1,"",tabA4[14,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Vuong Test of Levels vs. Changes:} \\\\")
latex_table_entry(file1,"\\hspace{.4cm}p-Value for null that parameter is in the overlapping region",tabA4[18,])
latex_table_entry(file1,"\\hspace{.4cm}p-Value for null that models fit equally well",tabA4[19,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & IV \\\\")
latex_table_entry(file1,"\\textbf{Number of Months}",tabA4[15,])
latex_table_entry(file1,"\\textbf{Number of Countries}",tabA4[16,])
latex_table_entry(file1,"\\textbf{R-Squared}",tabA4[17,])
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Elections Outcomes, the Economy, and Sentiment, Changes in Inflation --- Heteroskedasticity robust standard errors in parentheses. Column (5) corrects the results of column (4) for measurement error using instrumental variables. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{VoteShCoefTableInfCh}")
latex_add_line(file1,"\\end{table}\n\n\n")

# table A.5 #

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2)  \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA5[1,])
latex_table_entry(file1,"",tabA5[2,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA5[3,])
latex_table_entry(file1,"",tabA5[4,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA5[5,])
latex_table_entry(file1,"",tabA5[6,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA5[7,])
latex_table_entry(file1,"",tabA5[8,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA5[9,])
latex_table_entry(file1,"",tabA5[10,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA5[11,])
latex_table_entry(file1,"",tabA5[12,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA5[13,])
latex_table_entry(file1,"",tabA5[14,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA5[15,])
latex_table_entry(file1,"",tabA5[16,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Inflation (changes):}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA5[17,])
latex_table_entry(file1,"",tabA5[18,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA5[19,])
latex_table_entry(file1,"",tabA5[20,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA5[21,])
latex_table_entry(file1,"",tabA5[22,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA5[23,])
latex_table_entry(file1,"",tabA5[24,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & IV \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis, Changes in Inflation --- Standard errors calculated using the delta method in parentheses. Results in columns (2) and (3) are calculated based on columns (1) and (2) of Table \\ref{VoteShCoefTableInfCh}, respectively. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTableInfCh}")
latex_add_line(file1,"\\end{table}\n\n\n")

## table A.6 ##

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA6[1,])
latex_table_entry(file1,"",tabA6[2,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA6[3,])
latex_table_entry(file1,"",tabA6[4,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA6[5,])
latex_table_entry(file1,"",tabA6[6,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA6[7,])
latex_table_entry(file1,"",tabA6[8,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA6[9,])
latex_table_entry(file1,"",tabA6[10,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA6[11,])
latex_table_entry(file1,"",tabA6[12,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA6[13,])
latex_table_entry(file1,"",tabA6[14,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA6[15,])
latex_table_entry(file1,"",tabA6[16,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Inflation:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA6[17,])
latex_table_entry(file1,"",tabA6[18,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA6[19,])
latex_table_entry(file1,"",tabA6[20,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA6[21,])
latex_table_entry(file1,"",tabA6[22,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA6[23,])
latex_table_entry(file1,"",tabA6[24,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & IV \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis (clustered standard errors) --- Standard errors calculated using the delta method in parentheses. Results in columns (1), (2), and (3) are calculated based on columns (3), (4), and (5) of Table \\ref{VoteShCoefTable}, respectively. The first column reports mediation analysis for vote share, the second reports mediation analysis for vote intent, and the third corrects the results in column (2) for measurement error using instrumental variables. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTableClustered}")
latex_add_line(file1,"\\end{table}\n\n\n")

# table A.7 #
latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA7[1,])
latex_table_entry(file1,"",tabA7[2,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA7[3,])
latex_table_entry(file1,"",tabA7[4,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA7[5,])
latex_table_entry(file1,"",tabA7[6,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA7[7,])
latex_table_entry(file1,"",tabA7[8,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA7[9,])
latex_table_entry(file1,"",tabA7[10,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA7[11,])
latex_table_entry(file1,"",tabA7[12,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA7[13,])
latex_table_entry(file1,"",tabA7[14,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA7[15,])
latex_table_entry(file1,"",tabA7[16,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Inflation:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA7[17,])
latex_table_entry(file1,"",tabA7[18,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA7[19,])
latex_table_entry(file1,"",tabA7[20,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA7[21,])
latex_table_entry(file1,"",tabA7[22,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA7[23,])
latex_table_entry(file1,"",tabA7[24,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & IV \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis (panel Newey-West standard errors) --- Standard errors calculated using the delta method in parentheses. Results in columns (1), (2), and (3) are calculated based on columns (3), (4), and (5) of Table \\ref{VoteShCoefTable}, respectively. The first column reports mediation analysis for vote share, the second reports mediation analysis for vote intent, and the third corrects the results in column (2) for measurement error using instrumental variables. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTableNW}")
latex_add_line(file1,"\\end{table}\n\n\n")

# table A.8 #
latex_add_line(file1,"\\begin{landscape}")
latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c c c c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) & (4) & (5) & (6) & (7) & (8) \\\\")
latex_add_line(file1,"\\textbf{Dependent Variable:} & \\multicolumn{4}{c}{Vote Share} & \\multicolumn{4}{c}{Vote Intent} \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Independent Variables:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Constant",tabA8[1,])
latex_table_entry(file1,"",tabA8[2,])
latex_table_entry(file1,"\\hspace{.4cm}Growth (yearly)",tabA8[3,])
latex_table_entry(file1,"",tabA8[4,])
latex_table_entry(file1,"\\hspace{.4cm}Unem. (yearly)",tabA8[5,])
latex_table_entry(file1,"",tabA8[6,])
latex_table_entry(file1,"\\hspace{.4cm}Unem. Ch. (yearly)",tabA8[9,])
latex_table_entry(file1,"",tabA8[10,])
latex_table_entry(file1,"\\hspace{.4cm}Inf. (yearly)",tabA8[7,])
latex_table_entry(file1,"",tabA8[8,])
latex_table_entry(file1,"\\hspace{.4cm}Inf. Ch. (yearly)",tabA8[11,])
latex_table_entry(file1,"",tabA8[12,])
latex_add_line(file1,"\\\\")
latex_table_entry(file1,"\\textbf{Number of Months}",tabA8[13,])
latex_table_entry(file1,"\\textbf{Number of Countries}",tabA8[14,])
latex_table_entry(file1,"\\textbf{R-Squared}",dig3(as.numeric(tabA8[15,])))
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Elections Outcomes and the Economy --- Heteroskedasticity robust standard errors in parentheses. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{VoteShTotTable}")
latex_add_line(file1,"\\end{table}\n\n\n")
latex_add_line(file1,"\\end{landscape}")

## table A.9 ##

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3)  \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA9[1,])
latex_table_entry(file1,"",tabA9[2,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA9[3,])
latex_table_entry(file1,"",tabA9[4,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA9[5,])
latex_table_entry(file1,"",tabA9[6,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA9[7,])
latex_table_entry(file1,"",tabA9[8,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA9[9,])
latex_table_entry(file1,"",tabA9[10,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA9[11,])
latex_table_entry(file1,"",tabA9[12,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA9[13,])
latex_table_entry(file1,"",tabA9[14,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA9[15,])
latex_table_entry(file1,"",tabA9[16,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Inflation:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA9[17,])
latex_table_entry(file1,"",tabA9[18,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA9[19,])
latex_table_entry(file1,"",tabA9[20,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA9[21,])
latex_table_entry(file1,"",tabA9[22,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA9[23,])
latex_table_entry(file1,"",tabA9[24,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & IV \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis, LDV Specification --- Standard errors calculated using the delta method in parentheses. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTableLDV}")
latex_add_line(file1,"\\end{table}\n\n\n")

## table A.10 ##


latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA10[1,])
latex_table_entry(file1,"",tabA10[2,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA10[3,])
latex_table_entry(file1,"",tabA10[4,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA10[5,])
latex_table_entry(file1,"",tabA10[6,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA10[7,])
latex_table_entry(file1,"",tabA10[8,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA10[9,])
latex_table_entry(file1,"",tabA10[10,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA10[11,])
latex_table_entry(file1,"",tabA10[12,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA10[13,])
latex_table_entry(file1,"",tabA10[14,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA10[15,])
latex_table_entry(file1,"",tabA10[16,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & IV \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis (no inflation) --- Standard errors calculated using the delta method in parentheses. Results in columns (1), (2), and (3) are calculated based on columns (3), (4), and (5) of Table \\ref{VoteShCoefTable}, respectively. The first column reports mediation analysis for vote share, the second reports mediation analysis for vote intent, and the third corrects the results in column (2) for measurement error using instrumental variables. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTableNoInflation}")
latex_add_line(file1,"\\end{table}\n\n\n")

## table A.11 ##

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA11[1,])
latex_table_entry(file1,"",tabA11[2,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA11[3,])
latex_table_entry(file1,"",tabA11[4,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA11[5,])
latex_table_entry(file1,"",tabA11[6,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA11[7,])
latex_table_entry(file1,"",tabA11[8,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA11[9,])
latex_table_entry(file1,"",tabA11[10,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA11[11,])
latex_table_entry(file1,"",tabA11[12,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA11[13,])
latex_table_entry(file1,"",tabA11[14,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA11[15,])
latex_table_entry(file1,"",tabA11[16,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Inflation:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA11[17,])
latex_table_entry(file1,"",tabA11[18,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA11[19,])
latex_table_entry(file1,"",tabA11[20,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA11[21,])
latex_table_entry(file1,"",tabA11[22,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA11[23,])
latex_table_entry(file1,"",tabA11[24,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & IV \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis (trimmed inflation) --- Standard errors calculated using the delta method in parentheses. Results in columns (1), (2), and (3) are calculated based on columns (3), (4), and (5) of Table \\ref{VoteShCoefTable}, respectively. The first column reports mediation analysis for vote share, the second reports mediation analysis for vote intent, and the third corrects the results in column (2) for measurement error using instrumental variables. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTableTrimmed}")
latex_add_line(file1,"\\end{table}\n\n\n")

## table A.12 ##

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3)  \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA12[1,])
latex_table_entry(file1,"",tabA12[2,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA12[3,])
latex_table_entry(file1,"",tabA12[4,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA12[5,])
latex_table_entry(file1,"",tabA12[6,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA12[7,])
latex_table_entry(file1,"",tabA12[8,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA12[9,])
latex_table_entry(file1,"",tabA12[10,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA12[11,])
latex_table_entry(file1,"",tabA12[12,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA12[13,])
latex_table_entry(file1,"",tabA12[14,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA12[15,])
latex_table_entry(file1,"",tabA12[16,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Inflation:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA12[17,])
latex_table_entry(file1,"",tabA12[18,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA12[19,])
latex_table_entry(file1,"",tabA12[20,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA12[21,])
latex_table_entry(file1,"",tabA12[22,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA12[23,])
latex_table_entry(file1,"",tabA12[24,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & IV \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis, Additional Control Variables --- Standard errors calculated using the delta method in parentheses. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTableControls}")
latex_add_line(file1,"\\end{table}\n\n\n")

## table A.13 ##

latex_table_start(file1,3)
latex_add_line(file1," & (1) & (2) & (3) \\\\")
latex_add_line(file1,"\\textbf{Dependent Variable:}  & PM Party & Junior Party & Gov. Party \\\\")
latex_add_line(file1," & Vote & Vote & Vote \\\\")
latex_model(file1,c("Intercept","Growth","Unemployment","Inflation"),tabA13,eqlabels=F,stats=c("Number of Elections","Number of Countries","R-Squared"),kstat=c(0,0,3))
latex_table_end(file1,caption="Economic Vote, Prime Minister's Party vs. Junior Members --- Heteroskedasticity robust standard errors in parentheses. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.",label="JuniorTable")

## table A.14 ##

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) & (4) & (5) \\\\")
latex_add_line(file1,"\\textbf{Dependent Variable:} & Vote Share & Vote Intent & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Independent Variables:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Constant",tabA14[1,])
latex_table_entry(file1,"",tabA14[2,])
latex_table_entry(file1,"\\hspace{.4cm}Growth (yearly)",tabA14[3,])
latex_table_entry(file1,"",tabA14[4,])
latex_table_entry(file1,"\\hspace{.4cm}Unem. (yearly)",tabA14[5,])
latex_table_entry(file1,"",tabA14[6,])
latex_table_entry(file1,"\\hspace{.4cm}Inf. (yearly)",tabA14[7,])
latex_table_entry(file1,"",tabA14[8,])
latex_table_entry(file1,"\\hspace{.4cm}Growth Sentiment",tabA14[9,])
latex_table_entry(file1,"",tabA14[10,])
latex_table_entry(file1,"\\hspace{.4cm}Unem. Sentiment",tabA14[11,])
latex_table_entry(file1,"",tabA14[12,])
latex_table_entry(file1,"\\hspace{.4cm}Inf. Sentiment",tabA14[13,])
latex_table_entry(file1,"",tabA14[14,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & OLS & OLS & IV \\\\")
latex_table_entry(file1,"\\textbf{Number of Months}",tabA14[15,])
latex_table_entry(file1,"\\textbf{Number of Countries}",tabA14[16,])
latex_table_entry(file1,"\\textbf{R-Squared}",tabA14[17,])
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Elections Outcomes for Parties in Government, the Economy, and Sentiment --- Heteroskedasticity robust standard errors in parentheses. Column (5) corrects the results of column (4) for measurement error using instrumental variables. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{VoteShCoefTableGov}")
latex_add_line(file1,"\\end{table}\n\n\n")

## table A.15 ##

latex_add_line(file1,"\\begin{table}\\centering")
latex_add_line(file1,"\\begin{footnotesize}")
latex_add_line(file1,"\\begin{tabular}{ l c c c} \\hline \\hline")
latex_add_line(file1," & (1) & (2) & (3) \\\\")
latex_add_line(file1," \\textbf{Dependent Variable:}  & Vote Share & Vote Intent & Vote Intent \\\\")
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Growth:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA15[1,])
latex_table_entry(file1,"",tabA15[2,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA15[3,])
latex_table_entry(file1,"",tabA15[4,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA15[5,])
latex_table_entry(file1,"",tabA15[6,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA15[7,])
latex_table_entry(file1,"",tabA15[8,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Unemployment:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA15[9,])
latex_table_entry(file1,"",tabA15[10,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA15[11,])
latex_table_entry(file1,"",tabA15[12,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA15[13,])
latex_table_entry(file1,"",tabA15[14,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA15[15,])
latex_table_entry(file1,"",tabA15[16,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Inflation:}\\\\")
latex_table_entry(file1,"\\hspace{.4cm}Total Effect",tabA15[17,])
latex_table_entry(file1,"",tabA15[18,])
latex_table_entry(file1,"\\hspace{.4cm}Direct Effect",tabA15[19,])
latex_table_entry(file1,"",tabA15[20,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect Effect",tabA15[21,])
latex_table_entry(file1,"",tabA15[22,])
latex_table_entry(file1,"\\hspace{.4cm}Indirect / Total Ratio",tabA15[23,])
latex_table_entry(file1,"",tabA15[24,])
latex_add_line(file1,"\\\\")
latex_add_line(file1,"\\textbf{Estimator} & OLS & OLS & IV \\\\")
latex_add_line(file1,"\\hline \\hline")
latex_add_line(file1,"\\end{tabular}")
latex_add_line(file1,"\\end{footnotesize}")
latex_add_line(file1,"\\caption{Mediation Analysis, Vote Share of Parties in Government --- Standard errors calculated using the delta method in parentheses. The first column reports mediation analysis for vote share, the second reports mediation analysis for vote intent, and the third corrects the results in column (2) for measurement error using instrumental variables. In each case, the direct effect measures the effect of the economic variable controlling for sentiment and the indirect effect measures the effect of the economic variable as mediated by sentiment. $^+ p<.10, ^{*} p<.05,  ^{**} p<.01, ^{***} p<.001$.}")
latex_add_line(file1,"\\label{MediationTableGov}")
latex_add_line(file1,"\\end{table}\n\n\n")

latex_doc_end(file1)
